Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

PyomoNLP scaling factors on sub-blocks #3295

Merged
merged 19 commits into from
Aug 16, 2024
Merged
Show file tree
Hide file tree
Changes from 15 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 9 additions & 4 deletions pyomo/contrib/pynumero/interfaces/pyomo_grey_box_nlp.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,7 @@
)
from pyomo.contrib.pynumero.interfaces.external_grey_box import ExternalGreyBoxBlock
from pyomo.contrib.pynumero.interfaces.nlp_projections import ProjectedNLP
from pyomo.core.base.suffix import SuffixFinder


# Todo: make some of the numpy arrays not writable from __init__
Expand Down Expand Up @@ -227,12 +228,16 @@ def __init__(self, pyomo_model):
need_scaling = True

self._primals_scaling = np.ones(self.n_primals())
scaling_suffix = pyomo_model.component('scaling_factor')
scaling_suffix_finder = SuffixFinder('scaling_factor')
for i, v in enumerate(self._pyomo_model_var_datas):
v_scaling = scaling_suffix_finder.find(v)
if v_scaling is not None:
need_scaling = True
self._primals_scaling[i] = v_scaling
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is a minor edge case, but should we enforce that we never use a suffix that is "outside the scope" (above) the pyomo_model that was used to construct the PyomoNLP? Maybe SuffixFinder should support a context argument.

Copy link
Contributor

@Robbybp Robbybp Jul 30, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@blnicho @jsiirola @michaelbynum Can I get a second opinion here? The question is: what should happen when a PyomoNLP is constructed from a subblock that has relevant scaling factors on the parent block:

import pyomo.environ as pyo
from pyomo.contrib.pynumero.interfaces.pyomo_nlp import PyomoNLP

m = pyo.ConcreteModel()
m.b = pyo.Block()
m.b.x = pyo.Var([1, 2], initialize={1: 100, 2: 20})

# Components so we don't have an empty NLP
m.b.eq = pyo.Constraint(expr=m.b.x[1]*m.b.x[2] == 2000)
m.b.obj = pyo.Objective(expr=m.b.x[1]**2 + m.b.x[2]**2)

m.scaling_factor = pyo.Suffix(direction=pyo.Suffix.EXPORT)
m.scaling_factor[m.b.x[1]] = 1e-2
m.scaling_factor[m.b.x[2]] = 1e-1

nlp = PyomoNLP(m.b)
scaling = nlp.get_primals_scaling()
print(scaling)

Current behavior:

None

Proposed by this PR:

[0.01, 0.1 ]

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think that I would advocate the old behavior (although I can see value in both approaches!).

My rationale is the following: if we define a "model" as all "active components reachable through active blocks contained within the reference block," then we should exclude Suffixes outside of the subtree, as Suffix is an active component.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm inclined to agree.

# maintain backwards compatibility
scaling_suffix = self._pyomo_model.component('scaling_factor')
if scaling_suffix and scaling_suffix.ctype is pyo.Suffix:
need_scaling = True
for i, v in enumerate(self._pyomo_model_var_datas):
if v in scaling_suffix:
self._primals_scaling[i] = scaling_suffix[v]

self._constraints_scaling = BlockVector(len(nlps))
for i, nlp in enumerate(nlps):
Expand Down
54 changes: 36 additions & 18 deletions pyomo/contrib/pynumero/interfaces/pyomo_nlp.py
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@
from ..sparse.block_matrix import BlockMatrix
from pyomo.contrib.pynumero.interfaces.ampl_nlp import AslNLP
from pyomo.contrib.pynumero.interfaces.nlp import NLP
from pyomo.core.base.suffix import SuffixFinder
from .external_grey_box import ExternalGreyBoxBlock


Expand Down Expand Up @@ -298,34 +299,47 @@ def get_inequality_constraint_indices(self, constraints):
# overloaded from NLP
def get_obj_scaling(self):
obj = self.get_pyomo_objective()
val = SuffixFinder('scaling_factor').find(obj)
# maintain backwards compatibility
scaling_suffix = self._pyomo_model.component('scaling_factor')
if scaling_suffix and scaling_suffix.ctype is pyo.Suffix:
if obj in scaling_suffix:
return scaling_suffix[obj]
return 1.0
return None
return 1.0 if val is None else val
else:
return val

# overloaded from NLP
def get_primals_scaling(self):
scaling_suffix_finder = SuffixFinder('scaling_factor')
primals_scaling = np.ones(self.n_primals())
ret = None
for i, v in enumerate(self.get_pyomo_variables()):
val = scaling_suffix_finder.find(v)
if val is not None:
primals_scaling[i] = val
ret = primals_scaling
# maintain backwards compatibility
scaling_suffix = self._pyomo_model.component('scaling_factor')
if scaling_suffix and scaling_suffix.ctype is pyo.Suffix:
primals_scaling = np.ones(self.n_primals())
for i, v in enumerate(self.get_pyomo_variables()):
if v in scaling_suffix:
primals_scaling[i] = scaling_suffix[v]
return primals_scaling
return None
else:
return ret

# overloaded from NLP
def get_constraints_scaling(self):
scaling_suffix_finder = SuffixFinder('scaling_factor')
constraints_scaling = np.ones(self.n_constraints())
ret = None
for i, c in enumerate(self.get_pyomo_constraints()):
val = scaling_suffix_finder.find(c)
if val is not None:
constraints_scaling[i] = val
ret = constraints_scaling
# maintain backwards compatibility
scaling_suffix = self._pyomo_model.component('scaling_factor')
if scaling_suffix and scaling_suffix.ctype is pyo.Suffix:
constraints_scaling = np.ones(self.n_constraints())
for i, c in enumerate(self.get_pyomo_constraints()):
if c in scaling_suffix:
constraints_scaling[i] = scaling_suffix[c]
return constraints_scaling
return None
else:
return ret

def extract_subvector_grad_objective(self, pyomo_variables):
"""Compute the gradient of the objective and return the entries
Expand Down Expand Up @@ -607,12 +621,16 @@ def __init__(self, pyomo_model):
need_scaling = True

self._primals_scaling = np.ones(self.n_primals())
scaling_suffix = self._pyomo_nlp._pyomo_model.component('scaling_factor')
scaling_suffix_finder = SuffixFinder('scaling_factor')
for i, v in enumerate(self.get_pyomo_variables()):
v_scaling = scaling_suffix_finder.find(v)
if v_scaling is not None:
need_scaling = True
self._primals_scaling[i] = v_scaling
# maintain backwards compatibility
scaling_suffix = self._pyomo_model.component('scaling_factor')
if scaling_suffix and scaling_suffix.ctype is pyo.Suffix:
need_scaling = True
for i, v in enumerate(self.get_pyomo_variables()):
if v in scaling_suffix:
self._primals_scaling[i] = scaling_suffix[v]

self._constraints_scaling = []
pyomo_nlp_scaling = self._pyomo_nlp.get_constraints_scaling()
Expand Down
19 changes: 19 additions & 0 deletions pyomo/contrib/pynumero/interfaces/tests/test_nlp.py
Original file line number Diff line number Diff line change
Expand Up @@ -699,6 +699,25 @@ def test_indices_methods(self):
dense_hess = hess.todense()
self.assertTrue(np.array_equal(dense_hess, expected_hess))

def test_subblock_scaling(self):
m = pyo.ConcreteModel()
m.b = b = pyo.Block()
b.x = pyo.Var(bounds=(5e-17, 5e-16), initialize=1e-16)
b.scaling_factor = pyo.Suffix(direction=pyo.Suffix.EXPORT)
b.scaling_factor[b.x] = 1e16

b.c = pyo.Constraint(rule=b.x == 1e-16)
b.scaling_factor[b.c] = 1e16

b.o = pyo.Objective(expr=b.x)
b.scaling_factor[b.o] = 1e16

nlp = PyomoNLP(m)

assert nlp.get_obj_scaling() == 1e16
assert nlp.get_primals_scaling()[0] == 1e16
assert nlp.get_constraints_scaling()[0] == 1e16

def test_no_objective(self):
m = pyo.ConcreteModel()
m.x = pyo.Var()
Expand Down
Loading