Add new AutoScaler and CustomScalerBase classes #1429
Conversation
| jac, nlp = get_jacobian(blk, scaled=True) | ||
|
|
||
| if con_list is None: | ||
| con_list = nlp.get_pyomo_equality_constraints() |
There was a problem hiding this comment.
In my experience, get_jacobian is slow because forming the Pynumero NLP is slow. If we have flowsheet-level variables, we end up getting a Jacobian for the flowsheet and every subcomponent. I'm not sure how slow forming the NLP is going to be for PropertyBlock 43 in the MEA column, but you're still at least doubling the work you're doing.
I would suggest getting the Jacobian and NLP once at the beginning of constraints_by_jacobian_norm, then just referring to rows and columns of the overall object as you want to scale constraints.
There was a problem hiding this comment.
I have been thinking about there best way to do this, and did realize I was probably getting the Jacobian more times than necessary. I was not aware of the issue of having higher-level components in the problem causing the entire problem to be emitted (and given how properties work that will be an issue). I am also thinking there might need to be some checks to deal with deactivated components (and that it might be possible to deactivate unnecessary components before getting the Jacobian to reduce the problem size).
There was a problem hiding this comment.
I think I've fixed this so this method will only be called once for any call to the scaling method.
constraints_by_jacobian_norm (line 156) collects all the constraints we want to scale and the top level block, and then only calls the autoscaler once.
| # Use scipy to get all the norms | ||
| # Should be more efficient that iterating in Python | ||
| axis = ( | ||
| 1 # Could make this an argument to also support variable-based norm scaling |
There was a problem hiding this comment.
Based on your comment, should we consider adding another method for variable-based norm scaling?
There was a problem hiding this comment.
I do not know if this makes sense or not - comments would be welcome (@dallan-keylogic ).
There was a problem hiding this comment.
If we scale variables first based on order-of-magnitude, we can then scale constraints based on row/column norm. However, norm-based variable scaling probably won't be useful for a fully unscaled model, because the most common Jacobian entries will be O(1).
For a partially-scaled model, though, filling in variable scaling factors based on column norms might be useful.
idaes/core/scaling/autoscaling.py
Outdated
|
|
||
|
|
||
| @document_kwargs_from_configdict(CONFIG) | ||
| class AutoScaler: |
There was a problem hiding this comment.
It'd be nice to have a method like report_scaling that would give you vars and constraints in one column, scaling factor values in the 2nd column, and the scaled value in the third column.
There was a problem hiding this comment.
The first two you can get by doing suffix.display(). The third part - having scaled value is easy enough for Vars, but what is the "scaled value" for a Constraint? Scaled residual is easy to get, but not always that meaningful (as the scaled Jacobian is often more meaningful).
| assert jacobian_cond(methane, scaled=True) == pytest.approx(9191, abs=1) | ||
|
|
||
| # Check for optimal solution | ||
| assert check_optimal_termination(results) | ||
| assert len(extreme_jacobian_rows(methane, scaled=True)) == 0 | ||
| assert len(extreme_jacobian_rows(methane, scaled=True)) == 0 |
There was a problem hiding this comment.
Actually, adding to my other suggestion on thinking about a reporting method for Autoscaler, these would also be good to show in report output (e.g., condition number, extreme jac rows/columns or at least the number of each).
I realize that the DiagnosticsToolbox already has the capability to report on some of these things, but just thinking it'd be nice to be able to apply scaling via AutoScaler and subsequently check what you really did to model scaling with a report method.
There was a problem hiding this comment.
That is what the Diagnostics tools are for - putting those methods here would just be duplication of code (and the diagnostics also tell you a lot more). The eventual documentation will highlight that you really need to use diagnostics, scaling and initialization tools together to get the best results.
There was a problem hiding this comment.
Will the diagnostics toolbox reference these new tools?
There was a problem hiding this comment.
@agarciadiego I am not sure there - we probably need to think about how and where we draw the line between the two. Up until now, the Diagnostics Toolbox has been about finding the issues, and hasn't said much about how to fix them (as that is not so easy to explain concisely). We can definitely mention it in the docs, and maybe we could have the display methods mention tools to help fix issues (where appropriate).
…into autoscaling
| * large Jacobian Condition Numbers (>1e10), | ||
| * variables with very large or very small values, | ||
| * variables with values close to zero, | ||
| * variables with values close to a bound (conditionally), | ||
| * constraints with mismatched terms, | ||
| * constraints with potential cancellations, | ||
| * variables and constraints with extreme Jacobian norms, | ||
| * extreme entries in Jacobian matrix (conditionally). |
There was a problem hiding this comment.
These are all issues conditionally, although extreme Jacobian entries is more conditional than the others.
I believe "constraints with potential cancellations" aren't typically scaling issues but model formulation issues. "variables with values close to zero" is redundant with "variables with very large or very small values". "Variables with values close to bounds" typically end up being variables close to zero. I don't recall an example with a variable close to a nonzero bound that turned out to be bad scaling.
You could rank order them:
- Variables and constraints with extreme Jacobian column/row norms
- Variables with very large or very small values
- Gap
- Variables close to bounds
- Constraints with mismatched terms
- Extreme Jacobian Entries
As a side note, an instructive exercise for constraints with mismatched terms is to imagine setting the mismatched term to zero. If the Jacobian is nonsingular, it's probably fine. If it becomes singular, it's probably an issue.
There was a problem hiding this comment.
For this documentation I deliberately didn't go into too much detail. There is definitely room for a more detailed discussion of all these issues, but for the basic user I just wanted to give a list of things to look for.
docs/explanations/scaling_toolbox/diagnosing_scaling_issues.rst
Outdated
Show resolved
Hide resolved
| Args: | ||
| model: parent model of submodel | ||
| submodel: name of submodel to be scaled as str | ||
| submodel_scalers: user provided dict of Scalers to use for submodels |
There was a problem hiding this comment.
We could use the actual submodel objects if we did a ComponentMap. I assume you're aware of the option and elected for strings and dictionaries for some performance reason.
There was a problem hiding this comment.
I think the reason for this was that I was thinking we did not have the objects yet, but as this is a user-facing function (and not a set of defaults) that is not the case. I will update this too.
| sdict["block_datas"] = {} | ||
| for bd in blk_or_suffix.values(): | ||
| sdict["block_datas"][bd.name] = _collect_block_suffixes( | ||
| bd, descend_into=descend_into | ||
| ) |
There was a problem hiding this comment.
I was going to ask about what happens if multiple scaling suffixes for the same variable is encountered on different blocks, but it looks like the full block structure is recreated in this function.
There was a problem hiding this comment.
For all of these tools, we assume that the IDAES Controlled Suffix will be on the parent block of the component. Pyomo does support putting scaling factors elsewhere, but that would require the user to manually set those (and thus they are responsible for managing them).
| if not overwrite and component in sfx: | ||
| _log.debug( | ||
| f"Existing scaling factor for {component.name} found and overwrite=False. " | ||
| "Scaling factor unchanged." | ||
| ) |
There was a problem hiding this comment.
Is there a way to suppress debug messages like this for the scaling routines only? When I turn on Debug logging, it's typically to look at the IPOPT output of an initialization solve. I guess that involves passing the output level to the initialization object, so I imagine that the scaling routines take an outlevel argument in the same way.
There was a problem hiding this comment.
Not quite - you should be able to control the logging level by each subpackage in IDAES. I.e. you can set different levels for idaes and idaes.core.scaling.
| if not isinstance(blk, (Block, BlockData)): | ||
| raise TypeError( | ||
| "report_scaling_factors: blk must be an instance of a Pyomo Block." | ||
| ) |
There was a problem hiding this comment.
Are ConcreteModels instances of Blocks?
| stream.write(f"{n + ' ' * (maxname - len(n))}{TAB}{i}\n") | ||
|
|
||
|
|
||
| def get_nominal_value(component): |
There was a problem hiding this comment.
This method is technically correct, but I'm not sure how useful this concept of a "nominal value" is going to be. We'll see how it gets used.
There was a problem hiding this comment.
This is the underlying function that is used by the expression walker to get the nominal value of an expression/constraint. It is primarily intended to be used in the expression walker (or limited similar circumstances).
| scaled = jacobian_cond(methane, scaled=True) | ||
|
|
||
| assert scaled < unscaled | ||
| assert scaled == pytest.approx(6.96238e15, rel=1e-5) |
There was a problem hiding this comment.
Once this gets merged, I want to see what's causing the condition number to still be so badly behaved.
There was a problem hiding this comment.
My guess is the fact that the model is uninitialized, but I gave it good initial guesses for the actual final state - i.e. there is a mismatch between the current uninitialized state and the scaling factors it was given. If we went and solved the model, I suspect the value would be much lower.
docs/explanations/scaling_toolbox/diagnosing_scaling_issues.rst
Outdated
Show resolved
Hide resolved
| Args: | ||
| constraint: constraint to set scaling factor for | ||
| scheme: method to apply for determining constraint scaling | ||
| 'harmonic_mean': (default) sum(1/abs(nominal value)) |
There was a problem hiding this comment.
Not for this PR, but schemes could be explained more in the documentation
docs/explanations/scaling_toolbox/diagnosing_scaling_issues.rst
Outdated
Show resolved
Hide resolved
Summary/Motivation:
This PR adds a draft of the new
Scalerclasses along with some general utility functions for the new scaling interface. I have generally created new functions even if older ones existed to make a clean break from the older API and assist with backward compatibility.Draft Documentation: Initial Outline of Documentation for Scaling Tools.docx
Changes proposed in this PR:
Legal Acknowledgement
By contributing to this software project, I agree to the following terms and conditions for my contribution: