Coverage for gpkit/nomials/math.py: 97%

Shortcuts on this page

r m x   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

435 statements  

1"""Signomial, Posynomial, Monomial, Constraint, & MonoEQCOnstraint classes""" 

2from collections import defaultdict 

3import numpy as np 

4from .core import Nomial 

5from .array import NomialArray 

6from .. import units 

7from ..constraints import SingleEquationConstraint 

8from ..globals import SignomialsEnabled 

9from ..small_classes import Numbers 

10from ..small_classes import HashVector, EMPTY_HV 

11from ..varkey import VarKey 

12from ..small_scripts import mag 

13from ..exceptions import (InvalidGPConstraint, InvalidPosynomial, 

14 PrimalInfeasible, DimensionalityError) 

15from .map import NomialMap 

16from .substitution import parse_subs 

17 

18 

19class Signomial(Nomial): 

20 """A representation of a Signomial. 

21 

22 Arguments 

23 --------- 

24 exps: tuple of dicts 

25 Exponent dicts for each monomial term 

26 cs: tuple 

27 Coefficient values for each monomial term 

28 require_positive: bool 

29 If True and Signomials not enabled, c <= 0 will raise ValueError 

30 

31 Returns 

32 ------- 

33 Signomial 

34 Posynomial (if the input has only positive cs) 

35 Monomial (if the input has one term and only positive cs) 

36 """ 

37 _c = _exp = None # pylint: disable=invalid-name 

38 

39 __hash__ = Nomial.__hash__ 

40 

41 def __init__(self, hmap=None, cs=1, require_positive=True): # pylint: disable=too-many-statements,too-many-branches 

42 if not isinstance(hmap, NomialMap): 

43 if hasattr(hmap, "hmap"): 

44 hmap = hmap.hmap 

45 elif isinstance(hmap, Numbers): 

46 hmap_ = NomialMap([(EMPTY_HV, mag(hmap))]) 

47 hmap_.units_of_product(hmap) 

48 hmap = hmap_ 

49 elif isinstance(hmap, dict): 

50 exp = HashVector({VarKey(k): v for k, v in hmap.items() if v}) 

51 hmap = NomialMap({exp: mag(cs)}) 

52 hmap.units_of_product(cs) 

53 else: 

54 raise ValueError("Nomial construction accepts only NomialMaps," 

55 " objects with an .hmap attribute, numbers," 

56 " or *(exp dict of strings, number).") 

57 super().__init__(hmap) 

58 if self.any_nonpositive_cs: 

59 if require_positive and not SignomialsEnabled: 

60 raise InvalidPosynomial("each c must be positive.") 

61 self.__class__ = Signomial 

62 elif len(self.hmap) == 1: 

63 self.__class__ = Monomial 

64 else: 

65 self.__class__ = Posynomial 

66 

67 def diff(self, var): 

68 """Derivative of this with respect to a Variable 

69 

70 Arguments 

71 --------- 

72 var : Variable key 

73 Variable to take derivative with respect to 

74 

75 Returns 

76 ------- 

77 Signomial (or Posynomial or Monomial) 

78 """ 

79 var = var.key 

80 if var not in self.vks: 

81 diff = NomialMap({EMPTY_HV: 0.0}) 

82 diff.units = None 

83 else: 

84 diff = self.hmap.diff(var) 

85 return Signomial(diff, require_positive=False) 

86 

87 def posy_negy(self): 

88 """Get the positive and negative parts, both as Posynomials 

89 

90 Returns 

91 ------- 

92 Posynomial, Posynomial: 

93 p_pos and p_neg in (self = p_pos - p_neg) decomposition, 

94 """ 

95 py, ny = NomialMap(), NomialMap() 

96 py.units, ny.units = self.units, self.units 

97 for exp, c in self.hmap.items(): 

98 if c > 0: 

99 py[exp] = c 

100 elif c < 0: 

101 ny[exp] = -c # -c to keep it a posynomial 

102 return Posynomial(py) if py else 0, Posynomial(ny) if ny else 0 

103 

104 def mono_approximation(self, x0): 

105 """Monomial approximation about a point x0 

106 

107 Arguments 

108 --------- 

109 x0 (dict): 

110 point to monomialize about 

111 

112 Returns 

113 ------- 

114 Monomial (unless self(x0) < 0, in which case a Signomial is returned) 

115 """ 

116 x0, _, _ = parse_subs(self.vks, x0) # use only varkey keys 

117 psub = self.hmap.sub(x0, self.vks, parsedsubs=True) 

118 if EMPTY_HV not in psub or len(psub) > 1: 

119 raise ValueError("Variables %s remained after substituting x0=%s" 

120 " into %s" % (psub, x0, self)) 

121 c0, = psub.values() 

122 c, exp = c0, HashVector() 

123 for vk in self.vks: 

124 val = float(x0[vk]) 

125 diff, = self.hmap.diff(vk).sub(x0, self.vks, 

126 parsedsubs=True).values() 

127 e = val*diff/c0 

128 if e: 

129 exp[vk] = e 

130 try: 

131 c /= val**e 

132 except OverflowError: 

133 raise OverflowError( 

134 "While approximating the variable %s with a local value of" 

135 " %s, %s/(%s**%s) overflowed. Try reducing the variable's" 

136 " value by changing its unit prefix, or specify x0 values" 

137 " for any free variables it's multiplied or divided by in" 

138 " the posynomial %s whose expected value is far from 1." 

139 % (vk, val, c, val, e, self)) 

140 hmap = NomialMap({exp: c}) 

141 hmap.units = self.units 

142 return Monomial(hmap) 

143 

144 def sub(self, substitutions, require_positive=True): 

145 """Returns a nomial with substitued values. 

146 

147 Usage 

148 ----- 

149 3 == (x**2 + y).sub({'x': 1, y: 2}) 

150 3 == (x).gp.sub(x, 3) 

151 

152 Arguments 

153 --------- 

154 substitutions : dict or key 

155 Either a dictionary whose keys are strings, Variables, or VarKeys, 

156 and whose values are numbers, or a string, Variable or Varkey. 

157 val : number (optional) 

158 If the substitutions entry is a single key, val holds the value 

159 require_positive : boolean (optional, default is True) 

160 Controls whether the returned value can be a Signomial. 

161 

162 Returns 

163 ------- 

164 Returns substituted nomial. 

165 """ 

166 return Signomial(self.hmap.sub(substitutions, self.vks), 

167 require_positive=require_positive) 

168 

169 def __le__(self, other): 

170 if isinstance(other, (Numbers, Signomial)): 

171 return SignomialInequality(self, "<=", other) 

172 return NotImplemented 

173 

174 def __ge__(self, other): 

175 if isinstance(other, (Numbers, Signomial)): 

176 return SignomialInequality(self, ">=", other) 

177 return NotImplemented 

178 

179 def __add__(self, other, rev=False): 

180 other_hmap = getattr(other, "hmap", None) 

181 if isinstance(other, Numbers): 

182 if other == 0: 

183 return Signomial(self.hmap) 

184 other_hmap = NomialMap({EMPTY_HV: mag(other)}) 

185 other_hmap.units_of_product(other) 

186 if other_hmap: 

187 astorder = (self, other) 

188 if rev: 

189 astorder = tuple(reversed(astorder)) 

190 out = Signomial(self.hmap + other_hmap) 

191 out.ast = ("add", astorder) 

192 return out 

193 return NotImplemented 

194 

195 def __mul__(self, other, rev=False): 

196 astorder = (self, other) 

197 if rev: 

198 astorder = tuple(reversed(astorder)) 

199 if isinstance(other, np.ndarray): 

200 s = NomialArray(self) 

201 s.ast = self.ast 

202 return s*other 

203 if isinstance(other, Numbers): 

204 if not other: # other is zero 

205 return other 

206 hmap = mag(other)*self.hmap 

207 hmap.units_of_product(self.hmap.units, other) 

208 out = Signomial(hmap) 

209 out.ast = ("mul", astorder) 

210 return out 

211 if isinstance(other, Signomial): 

212 hmap = NomialMap() 

213 for exp_s, c_s in self.hmap.items(): 

214 for exp_o, c_o in other.hmap.items(): 

215 exp = exp_s + exp_o 

216 new, accumulated = c_s*c_o, hmap.get(exp, 0) 

217 if new != -accumulated: 

218 hmap[exp] = accumulated + new 

219 elif accumulated: 

220 del hmap[exp] 

221 hmap.units_of_product(self.hmap.units, other.hmap.units) 

222 out = Signomial(hmap) 

223 out.ast = ("mul", astorder) 

224 return out 

225 return NotImplemented 

226 

227 def __truediv__(self, other): 

228 "Support the / operator in Python 2.x" 

229 if isinstance(other, Numbers): 

230 out = self*other**-1 

231 out.ast = ("div", (self, other)) 

232 return out 

233 if isinstance(other, Monomial): 

234 return other.__rtruediv__(self) 

235 return NotImplemented 

236 

237 def __pow__(self, expo): 

238 if isinstance(expo, int) and expo >= 0: 

239 p = 1 

240 while expo > 0: 

241 p *= self 

242 expo -= 1 

243 p.ast = ("pow", (self, expo)) 

244 return p 

245 return NotImplemented 

246 

247 def __neg__(self): 

248 if SignomialsEnabled: # pylint: disable=using-constant-test 

249 out = -1*self 

250 out.ast = ("neg", self) 

251 return out 

252 return NotImplemented 

253 

254 def __sub__(self, other): 

255 return self + -other if SignomialsEnabled else NotImplemented # pylint: disable=using-constant-test 

256 

257 def __rsub__(self, other): 

258 return other + -self if SignomialsEnabled else NotImplemented # pylint: disable=using-constant-test 

259 

260 def chop(self): 

261 "Returns a list of monomials in the signomial." 

262 monmaps = [NomialMap({exp: c}) for exp, c in self.hmap.items()] 

263 for monmap in monmaps: 

264 monmap.units = self.hmap.units 

265 return [Monomial(monmap) for monmap in sorted(monmaps, key=str)] 

266 

267class Posynomial(Signomial): 

268 "A Signomial with strictly positive cs" 

269 

270 __hash__ = Signomial.__hash__ 

271 

272 def __le__(self, other): 

273 if isinstance(other, Numbers + (Monomial,)): 

274 return PosynomialInequality(self, "<=", other) 

275 return NotImplemented 

276 

277 # Posynomial.__ge__ falls back on Signomial.__ge__ 

278 

279 def mono_lower_bound(self, x0): 

280 """Monomial lower bound at a point x0 

281 

282 Arguments 

283 --------- 

284 x0 (dict): 

285 point to make lower bound exact 

286 

287 Returns 

288 ------- 

289 Monomial 

290 """ 

291 return self.mono_approximation(x0) 

292 

293 

294class Monomial(Posynomial): 

295 "A Posynomial with only one term" 

296 

297 __hash__ = Posynomial.__hash__ 

298 

299 @property 

300 def exp(self): 

301 "Creates exp or returns a cached exp" 

302 if not self._exp: 

303 self._exp, = self.hmap.keys() # pylint: disable=attribute-defined-outside-init 

304 return self._exp 

305 

306 @property 

307 def c(self): # pylint: disable=invalid-name 

308 "Creates c or returns a cached c" 

309 if not self._c: 

310 self._c, = self.cs # pylint: disable=attribute-defined-outside-init, invalid-name 

311 return self._c 

312 

313 def __rtruediv__(self, other): 

314 "Divide other by this Monomial" 

315 if isinstance(other, Numbers + (Signomial,)): 

316 out = other * self**-1 

317 out.ast = ("div", (other, self)) 

318 return out 

319 return NotImplemented 

320 

321 def __pow__(self, expo): 

322 if isinstance(expo, Numbers): 

323 (exp, c), = self.hmap.items() 

324 exp = exp*expo if expo else EMPTY_HV 

325 hmap = NomialMap({exp: c**expo}) 

326 if expo and self.hmap.units: 

327 hmap.units = self.hmap.units**expo 

328 else: 

329 hmap.units = None 

330 out = Monomial(hmap) 

331 out.ast = ("pow", (self, expo)) 

332 return out 

333 return NotImplemented 

334 

335 def __eq__(self, other): 

336 if isinstance(other, MONS): 

337 try: # if both are monomials, return a constraint 

338 return MonomialEquality(self, other) 

339 except (DimensionalityError, ValueError) as e: 

340 print("Infeasible monomial equality: %s" % e) 

341 return False 

342 return super().__eq__(other) 

343 

344 def __ge__(self, other): 

345 if isinstance(other, Numbers + (Posynomial,)): 

346 return PosynomialInequality(self, ">=", other) 

347 # elif isinstance(other, np.ndarray): 

348 # return other.__le__(self, rev=True) 

349 return NotImplemented 

350 

351 # Monomial.__le__ falls back on Posynomial.__le__ 

352 

353 def mono_approximation(self, x0): 

354 return self 

355 

356 

357MONS = Numbers + (Monomial,) 

358 

359 

360####################################################### 

361####### CONSTRAINTS ################################### 

362####################################################### 

363 

364 

365class ScalarSingleEquationConstraint(SingleEquationConstraint): 

366 "A SingleEquationConstraint with scalar left and right sides." 

367 generated_by = v_ss = parent = None 

368 bounded = meq_bounded = {} 

369 

370 def __init__(self, left, oper, right): 

371 lr = [left, right] 

372 self.vks = set() 

373 for i, sig in enumerate(lr): 

374 if isinstance(sig, Signomial): 

375 for exp in sig.hmap: 

376 self.vks.update(exp) 

377 else: 

378 lr[i] = Signomial(sig) 

379 from .. import NamedVariables 

380 self.lineage = tuple(NamedVariables.lineage) 

381 super().__init__(lr[0], oper, lr[1]) 

382 

383 def relaxed(self, relaxvar): 

384 "Returns the relaxation of the constraint in a list." 

385 if self.oper == ">=": 

386 return [relaxvar*self.left >= self.right] 

387 if self.oper == "<=": 

388 return [self.left <= relaxvar*self.right] 

389 if self.oper == "=": 

390 return [self.left <= relaxvar*self.right, 

391 relaxvar*self.left >= self.right] 

392 raise ValueError( 

393 "Constraint %s had unknown operator %s." % self.oper, self) 

394 

395 

396# pylint: disable=too-many-instance-attributes, invalid-unary-operand-type 

397class PosynomialInequality(ScalarSingleEquationConstraint): 

398 """A constraint of the general form monomial >= posynomial 

399 Stored in the posylt1_rep attribute as a single Posynomial (self <= 1) 

400 Usually initialized via operator overloading, e.g. cc = (y**2 >= 1 + x) 

401 """ 

402 feastol = 1e-3 

403 # NOTE: follows .check_result's max default, but 1e-3 seems a bit lax... 

404 

405 def __init__(self, left, oper, right): 

406 ScalarSingleEquationConstraint.__init__(self, left, oper, right) 

407 if self.oper == "<=": 

408 p_lt, m_gt = self.left, self.right 

409 elif self.oper == ">=": 

410 m_gt, p_lt = self.left, self.right 

411 else: 

412 raise ValueError("operator %s is not supported." % self.oper) 

413 

414 self.unsubbed = self._gen_unsubbed(p_lt, m_gt) 

415 self.bounded = set() 

416 for p in self.unsubbed: 

417 for exp in p.hmap: 

418 for vk, x in exp.items(): 

419 self.bounded.add((vk, "upper" if x > 0 else "lower")) 

420 

421 def _simplify_posy_ineq(self, hmap, pmap=None, fixed=None): 

422 "Simplify a posy <= 1 by moving constants to the right side." 

423 if EMPTY_HV not in hmap: 

424 return hmap 

425 coeff = 1 - hmap[EMPTY_HV] 

426 if pmap is not None: # note constant term's mmap 

427 const_idx = list(hmap.keys()).index(EMPTY_HV) 

428 self.const_mmap = self.pmap.pop(const_idx) # pylint: disable=attribute-defined-outside-init 

429 self.const_coeff = coeff # pylint: disable=attribute-defined-outside-init 

430 if coeff >= -self.feastol and len(hmap) == 1: 

431 return None # a tautological monomial! 

432 if coeff < -self.feastol: 

433 msg = "'%s' is infeasible by %.2g%%" % (self, -coeff*100) 

434 if fixed: 

435 msg += " after substituting %s." % fixed 

436 raise PrimalInfeasible(msg) 

437 scaled = hmap/coeff 

438 scaled.units = hmap.units 

439 del scaled[EMPTY_HV] 

440 return scaled 

441 

442 def _gen_unsubbed(self, p_lt, m_gt): 

443 """Returns the unsubstituted posys <= 1. 

444 

445 Parameters 

446 ---------- 

447 p_lt : posynomial 

448 the left-hand side of (posynomial < monomial) 

449 

450 m_gt : monomial 

451 the right-hand side of (posynomial < monomial) 

452 

453 """ 

454 try: 

455 m_exp, = m_gt.hmap.keys() 

456 m_c, = m_gt.hmap.values() 

457 except ValueError: 

458 raise TypeError("greater-than side '%s' is not monomial." % m_gt) 

459 m_c *= units.of_division(m_gt, p_lt) 

460 hmap = p_lt.hmap.copy() 

461 for exp in list(hmap): 

462 hmap[exp-m_exp] = hmap.pop(exp)/m_c 

463 hmap = self._simplify_posy_ineq(hmap) 

464 return [Posynomial(hmap)] if hmap else [] 

465 

466 def as_hmapslt1(self, substitutions): 

467 "Returns the posys <= 1 representation of this constraint." 

468 out = [] 

469 for posy in self.unsubbed: 

470 fixed, _, _ = parse_subs(posy.vks, substitutions, clean=True) 

471 hmap = posy.hmap.sub(fixed, posy.vks, parsedsubs=True) 

472 self.pmap = hmap.mmap(posy.hmap) # pylint: disable=attribute-defined-outside-init 

473 del hmap.expmap, hmap.csmap # needed only for the mmap call above 

474 hmap = self._simplify_posy_ineq(hmap, self.pmap, fixed) 

475 if hmap is not None: 

476 if any(c <= 0 for c in hmap.values()): 

477 raise InvalidGPConstraint("'%s' became Signomial after sub" 

478 "stituting %s" % (self, fixed)) 

479 hmap.parent = self 

480 out.append(hmap) 

481 return out 

482 

483 def sens_from_dual(self, la, nu, _): 

484 "Returns the variable/constraint sensitivities from lambda/nu" 

485 presub, = self.unsubbed 

486 if hasattr(self, "pmap"): 

487 nu_ = np.zeros(len(presub.hmap)) 

488 for i, mmap in enumerate(self.pmap): 

489 for idx, percentage in mmap.items(): 

490 nu_[idx] += percentage*nu[i] 

491 del self.pmap # not needed after dual has been derived 

492 if hasattr(self, "const_mmap"): 

493 scale = (1-self.const_coeff)/self.const_coeff 

494 for idx, percentage in self.const_mmap.items(): 

495 nu_[idx] += percentage * la*scale 

496 del self.const_mmap # not needed after dual has been derived 

497 nu = nu_ 

498 self.v_ss = HashVector() 

499 if self.parent: 

500 self.parent.v_ss = self.v_ss 

501 if self.generated_by: 

502 self.generated_by.v_ss = self.v_ss 

503 for nu_i, exp in zip(nu, presub.hmap): 

504 for vk, x in exp.items(): 

505 self.v_ss[vk] = nu_i*x + self.v_ss.get(vk, 0) 

506 return self.v_ss, la 

507 

508 

509class MonomialEquality(PosynomialInequality): 

510 "A Constraint of the form Monomial == Monomial." 

511 oper = "=" 

512 

513 def __init__(self, left, right): 

514 # pylint: disable=super-init-not-called,non-parent-init-called 

515 ScalarSingleEquationConstraint.__init__(self, left, self.oper, right) 

516 self.unsubbed = self._gen_unsubbed(self.left, self.right) 

517 self.bounded = set() 

518 self.meq_bounded = {} 

519 self._las = [] 

520 if self.unsubbed and len(self.vks) > 1: 

521 exp, = self.unsubbed[0].hmap 

522 for key, e in exp.items(): 

523 s_e = np.sign(e) 

524 ubs = frozenset((k, "upper" if np.sign(e) != s_e else "lower") 

525 for k, e in exp.items() if k != key) 

526 lbs = frozenset((k, "lower" if np.sign(e) != s_e else "upper") 

527 for k, e in exp.items() if k != key) 

528 self.meq_bounded[(key, "upper")] = frozenset([ubs]) 

529 self.meq_bounded[(key, "lower")] = frozenset([lbs]) 

530 

531 def _gen_unsubbed(self, left, right): # pylint: disable=arguments-differ 

532 "Returns the unsubstituted posys <= 1." 

533 unsubbed = PosynomialInequality._gen_unsubbed 

534 l_over_r = unsubbed(self, left, right) 

535 r_over_l = unsubbed(self, right, left) 

536 return l_over_r + r_over_l 

537 

538 def as_hmapslt1(self, substitutions): 

539 "Tags posynomials for dual feasibility checking" 

540 out = super().as_hmapslt1(substitutions) 

541 for h in out: 

542 h.from_meq = True # pylint: disable=attribute-defined-outside-init 

543 return out 

544 

545 def __bool__(self): 

546 'A constraint not guaranteed to be satisfied evaluates as "False".' 

547 return bool(self.left.c == self.right.c 

548 and self.left.exp == self.right.exp) 

549 

550 def sens_from_dual(self, la, nu, _): 

551 "Returns the variable/constraint sensitivities from lambda/nu" 

552 self._las.append(la) 

553 if len(self._las) == 1: 

554 return {}, 0 

555 la = self._las[0] - self._las[1] 

556 self._las = [] 

557 exp, = self.unsubbed[0].hmap 

558 self.v_ss = exp*la 

559 return self.v_ss, la 

560 

561 

562class SignomialInequality(ScalarSingleEquationConstraint): 

563 """A constraint of the general form posynomial >= posynomial 

564 

565 Stored at .unsubbed[0] as a single Signomial (0 >= self)""" 

566 

567 def __init__(self, left, oper, right): 

568 ScalarSingleEquationConstraint.__init__(self, left, oper, right) 

569 if not SignomialsEnabled: 

570 raise TypeError("Cannot initialize SignomialInequality" 

571 " outside of a SignomialsEnabled environment.") 

572 if self.oper == "<=": 

573 plt, pgt = self.left, self.right 

574 elif self.oper == ">=": 

575 pgt, plt = self.left, self.right 

576 else: 

577 raise ValueError("operator %s is not supported." % self.oper) 

578 self.unsubbed = [plt - pgt] 

579 self.bounded = self.as_gpconstr({}).bounded 

580 

581 def as_hmapslt1(self, substitutions): 

582 "Returns the posys <= 1 representation of this constraint." 

583 siglt0, = self.unsubbed 

584 siglt0 = siglt0.sub(substitutions, require_positive=False) 

585 posy, negy = siglt0.posy_negy() 

586 if posy is 0: # pylint: disable=literal-comparison 

587 print("Warning: SignomialConstraint %s became the tautological" 

588 " constraint 0 <= %s after substitution." % (self, negy)) 

589 return [] 

590 if negy is 0: # pylint: disable=literal-comparison 

591 raise ValueError("%s became the infeasible constraint %s <= 0" 

592 " after substitution." % (self, posy)) 

593 if hasattr(negy, "cs") and len(negy.cs) > 1: 

594 raise InvalidGPConstraint( 

595 "%s did not simplify to a PosynomialInequality; try calling" 

596 " `.localsolve` instead of `.solve` to form your Model as a" 

597 " SequentialGeometricProgram." % self) 

598 # all but one of the negy terms becomes compatible with the posy 

599 p_ineq = PosynomialInequality(posy, "<=", negy) 

600 p_ineq.parent = self 

601 siglt0_us, = self.unsubbed 

602 siglt0_hmap = siglt0_us.hmap.sub(substitutions, siglt0_us.vks) 

603 negy_hmap = NomialMap() 

604 posy_hmaps = defaultdict(NomialMap) 

605 for o_exp, exp in siglt0_hmap.expmap.items(): 

606 if exp == negy.exp: 

607 negy_hmap[o_exp] = -siglt0_us.hmap[o_exp] 

608 else: 

609 posy_hmaps[exp-negy.exp][o_exp] = siglt0_us.hmap[o_exp] 

610 # pylint: disable=attribute-defined-outside-init 

611 self._mons = [Monomial(NomialMap({k: v})) 

612 for k, v in (posy/negy).hmap.items()] 

613 self._negysig = Signomial(negy_hmap, require_positive=False) 

614 self._coeffsigs = {exp: Signomial(hmap, require_positive=False) 

615 for exp, hmap in posy_hmaps.items()} 

616 self._sigvars = {exp: (list(self._negysig.vks) 

617 + list(sig.vks)) 

618 for exp, sig in self._coeffsigs.items()} 

619 return p_ineq.as_hmapslt1(substitutions) 

620 

621 def sens_from_dual(self, la, nu, result): 

622 """ We want to do the following chain: 

623 dlog(Obj)/dlog(monomial[i]) = nu[i] 

624 * dlog(monomial)/d(monomial) = 1/(monomial value) 

625 * d(monomial)/d(var) = see below 

626 * d(var)/dlog(var) = var 

627 = dlog(Obj)/dlog(var) 

628 each final monomial is really 

629 (coeff signomial)/(negy signomial) 

630 and by the chain rule d(monomial)/d(var) = 

631 d(coeff)/d(var)*1/negy + d(1/negy)/d(var)*coeff 

632 = d(coeff)/d(var)*1/negy - d(negy)/d(var)*coeff*1/negy**2 

633 """ 

634 # pylint: disable=too-many-locals, attribute-defined-outside-init 

635 

636 # pylint: disable=no-member 

637 def subval(posy): 

638 "Substitute solution into a posynomial and return the result" 

639 hmap = posy.sub(result["variables"], 

640 require_positive=False).hmap 

641 (key, value), = hmap.items() 

642 assert not key # constant 

643 return value 

644 

645 self.v_ss = {} 

646 invnegy_val = 1/subval(self._negysig) 

647 for i, nu_i in enumerate(nu): 

648 mon = self._mons[i] 

649 inv_mon_val = 1/subval(mon) 

650 coeff = self._coeffsigs[mon.exp] 

651 for var in self._sigvars[mon.exp]: 

652 d_mon_d_var = (subval(coeff.diff(var))*invnegy_val 

653 - (subval(self._negysig.diff(var)) 

654 * subval(coeff) * invnegy_val**2)) 

655 var_val = result["variables"][var] 

656 sens = (nu_i*inv_mon_val*d_mon_d_var*var_val) 

657 assert isinstance(sens, float) 

658 self.v_ss[var] = sens + self.v_ss.get(var, 0) 

659 return self.v_ss, la 

660 

661 def as_gpconstr(self, x0): 

662 "Returns GP-compatible approximation at x0" 

663 siglt0, = self.unsubbed 

664 posy, negy = siglt0.posy_negy() 

665 # default guess of 1.0 for unspecified negy variables 

666 x0 = {vk: x0.get(vk, 1) for vk in negy.vks} 

667 pconstr = PosynomialInequality(posy, "<=", negy.mono_lower_bound(x0)) 

668 pconstr.generated_by = self 

669 return pconstr 

670 

671 

672class SingleSignomialEquality(SignomialInequality): 

673 "A constraint of the general form posynomial == posynomial" 

674 

675 def __init__(self, left, right): 

676 SignomialInequality.__init__(self, left, "<=", right) 

677 self.oper = "=" 

678 self.meq_bounded = self.as_gpconstr({}).meq_bounded 

679 

680 def as_hmapslt1(self, substitutions): 

681 "SignomialEquality is never considered GP-compatible" 

682 raise InvalidGPConstraint(self) 

683 

684 def as_gpconstr(self, x0): 

685 "Returns GP-compatible approximation at x0" 

686 siglt0, = self.unsubbed 

687 posy, negy = siglt0.posy_negy() 

688 # default guess of 1.0 for unspecified negy variables 

689 x0 = {vk: x0.get(vk, 1) for vk in siglt0.vks} 

690 mec = (posy.mono_lower_bound(x0) == negy.mono_lower_bound(x0)) 

691 mec.generated_by = self 

692 return mec