Coverage for backend/ahuora-builder/src/ahuora_builder/flowsheet_manager.py: 85%
345 statements
« prev ^ index » next coverage.py v7.10.7, created at 2026-07-22 05:22 +0000
1from collections import defaultdict, deque
2import re
3from pyomo.network import SequentialDecomposition, Port
4from pyomo.environ import (
5 Component,
6 ConcreteModel,
7 TransformationFactory,
8 SolverFactory,
9 Block,
10 Expression,
11 Constraint,
12 Objective,
13 minimize,
14 assert_optimal_termination,
15 ScalarVar,
16 units as pyunits,
17 TerminationCondition
18)
19from pyomo.environ import assert_optimal_termination, check_optimal_termination
20from pyomo.core.base.constraint import ScalarConstraint, IndexedConstraint
21from idaes.core import FlowsheetBlock
22from .model_statistics_optimizations import (
23 degrees_of_freedom,
24 install_model_statistics_optimizations,
25)
27install_model_statistics_optimizations()
29from idaes.core.util.model_statistics import report_statistics
30from idaes.core.util.exceptions import InitializationError
31import idaes.logger as idaeslog
32from ahuora_builder_types import FlowsheetSchema
33from ahuora_builder_types.flowsheet_schema import SolvedFlowsheetSchema, UnitDiagnosticsFinding
34from .property_package_manager import PropertyPackageManager
35from .port_manager import PortManager
36from .arc_manager import ArcManager
37from .tear_manager import TearManager
38from .unit_model_manager import UnitModelManager
39from .methods.adapter_library import AdapterLibrary
40from .methods.adapter import (
41 serialize_properties_map,
42 deactivate_fixed_guesses,
43 deactivate_component,
44 deactivate_components,
45 add_corresponding_constraint,
46)
47from .methods.expression_parsing import parse_expression, ExpressionParsingError
48from .timing import start_timing
49from .methods.units_handler import get_value, get_attached_unit, check_units_equivalent
50from .methods.scaling_suffix import sanitize_scaling_suffix
51from .properties_manager import PropertiesManager
52from .custom.energy.power_property_package import PowerParameterBlock
53from .custom.energy.ac_property_package import acParameterBlock
54from .custom.energy.transformer_property_package import transformerParameterBlock
55from .custom.watertap.reverse_osmosis_0d import ReverseOsmosis0DData
56from pyomo.core.base.units_container import units, _PyomoUnit
57from idaes.core.util.model_serializer import StoreSpec, from_json, to_json
58from .diagnostics.infeasibilities import get_top_infeasibilities
60# CSTR Imports from example flowsheet.
61# To be used as placeholders until bespoke functions can be developed.
62from .custom import hda_reaction as reaction_props
63from .custom.hda_ideal_VLE import HDAParameterBlock
64from ahuora_builder.properties_manager import PropertyComponent
65from ahuora_property_packages.build_package import build_package
66from .methods.change_detection import get_current_values, detect_changes, print_changes
67from .methods.pyomo_expressions import dimensionless_objective_expression, sum_pyomo_expressions
69_log = idaeslog.getLogger(__name__)
71WATERTAP_SOLVER_WARNING = (
72 "Reverse Osmosis (0D) is a WaterTAP unit and should be solved with the "
73 "Ipopt Watertap scenario solver. Plain Ipopt can converge to a locally "
74 "infeasible point for this model even when the inputs are valid."
75)
78def _get_watertap_solver_diagnostics(
79 *,
80 solver_option: str | None,
81 unit_models,
82 properties_map: PropertiesManager,
83) -> list[UnitDiagnosticsFinding]:
84 """
85 Warn when a unit has a known dependency on the WaterTAP Ipopt wrapper.
87 This is intentionally solver-specific rather than part of the unit's
88 ``diagnose`` method, because unit diagnostics do not know which scenario
89 solver the user selected.
90 """
91 if solver_option != "ipopt":
92 return []
94 diagnostics: list[UnitDiagnosticsFinding] = []
95 for unit_model in unit_models:
96 if isinstance(unit_model, ReverseOsmosis0DData):
97 diagnostics.append(
98 UnitDiagnosticsFinding(
99 severity="warning",
100 property_value_id=properties_map.get_id_by_component(
101 unit_model.A_comp
102 ),
103 message=WATERTAP_SOLVER_WARNING,
104 )
105 )
106 return diagnostics
108# from amplpy import modules
109# Import required to allow the library to set the PATH and allow conopt to be found.
112def build_flowsheet(dynamic=False,time_set=[0]):
113 """
114 Builds a flowsheet block
115 """
116 # create the model and the flowsheet block
117 model = ConcreteModel()
118 model.fs = FlowsheetBlock(dynamic=dynamic, time_set=time_set, time_units=units.s)
119 model.fs.guess_vars = []
120 model.fs.controlled_vars = Block()
121 # Always add the power property package (as there's only one and there's no different property package types)
122 # TODO: should this be on-demand?
123 model.fs.power_pp = PowerParameterBlock()
124 model.fs.ac_pp = acParameterBlock()
125 model.fs.tr_pp = transformerParameterBlock()
126 # properties map: { id: pyomo variable or expression }
127 # used to map symbols for a sympy expression
128 model.fs.properties_map = PropertiesManager()
129 # list of component constraints to add later
130 model.fs.constraint_exprs = []
132 # Placholder property packages for CSTR reactor.
133 model.fs.BTHM_params = HDAParameterBlock()
134 #Hard-coded peng-robinson package for reactor:
135 model.fs.peng_robinson = build_package("peng-robinson",["benzene", "toluene", "hydrogen", "methane"], ["Liq","Vap"])
137 # Reaction package for the HDA reaction
138 model.fs.reaction_params = reaction_props.HDAReactionParameterBlock(
139 property_package=model.fs.peng_robinson
140 )
142 return model
145class FlowsheetManager:
146 """
147 Manages the flowsheet, including the property packages, unit models, ports, arcs, and tears
148 Includes methods to load, initialise, and solve the flowsheet
149 """
151 def __init__(self, schema: FlowsheetSchema) -> None:
152 """
153 Stores all relevant information about the flowsheet, without actually loading it
154 """
155 self.timing = start_timing()
156 self.timing.add_timing("initialise_flowsheet_manager")
158 self.model = build_flowsheet(dynamic=schema.dynamic, time_set=schema.time_set)
160 self.schema = schema
161 # Add property packages first, so that unit models can use them
162 self.property_packages = PropertyPackageManager(self)
163 # Add the port manager, so the unit models can register their ports
164 self.ports = PortManager()
165 # Add unit models
166 self.unit_models = UnitModelManager(self)
167 # Add arcs to connect the unit models together
168 self.arcs = ArcManager(self)
169 # set certain arcs as tears
170 self.tears = TearManager(self)
171 self.properties_map: PropertiesManager
173 def load(self) -> None:
174 """
175 Parses the schema and loads the model
176 """
177 self.timing.step_into("load_flowsheet")
178 # Load property packages first, so that unit models can use them
179 self.timing.add_timing("load_property_packages")
180 self.property_packages.load()
181 self.timing.step_into("load_unit_models")
182 self.unit_models.load()
183 self.timing.step_out()
184 # no need to load ports seperately, they are loaded by the unit models
185 # Load arcs to connect the unit models together
186 self.arcs.load()
187 # load any expressions
188 self.load_specs()
189 # if dynamics, apply finite difference transformation.
190 if self.schema.dynamic:
191 print("performing finite difference with", len(self.model.fs.time), "time steps")
192 TransformationFactory("dae.finite_difference").apply_to(
193 self.model.fs,
194 nfe=len(self.model.fs.time)-1, # Number of finite elements to use for discretization. We aren't adding any extra steps as our constraints dont work for that.
195 wrt=self.model.fs.time,
196 scheme="BACKWARD"
197 )
199 self.properties_map = self.model.fs.properties_map
200 self.timing.step_out()
202 def load_specs(self) -> None:
203 """
204 Loads expressions from the schema
205 """
206 fs = self.model.fs
207 specs = self.schema.expressions or []
209 ## Build dependency tree and sort specs before loading.
210 dependencies, result_ids = self.build_spec_dependency_tree(specs)
211 sorted_result_ids = self.topological_sort(dependencies, result_ids)
213 # load the specs (expressions within specifications tab)
214 for result_id in sorted_result_ids:
215 spec_config = next(spec for spec in specs if spec["id"] == result_id)
216 expression_str = spec_config["expression"]
217 try:
218 component_name = f"{spec_config['name']}_{spec_config['id']}"
219 def expression_rule(blk, time_index):
220 return parse_expression(expression_str, fs,time_index)
222 component = Expression(fs.time, rule=expression_rule)
223 fs.add_component(component_name, component)
224 fs.properties_map.add(
225 spec_config["id"], component, component.name, unknown_units=True
226 )
227 except ExpressionParsingError as e:
228 raise ExpressionParsingError(f"{e} when parsing expression '{expression_str}' for {component_name}: ")
230 # load constraints (expressions for specific property infos)
231 # can only handle equality constraints for now
232 for component, expr_str, id in fs.constraint_exprs:
233 # get the time index
234 #for time_index in fs.time:
235 if component.index_set().dimen != 0 and component.index_set() != fs.time: 235 ↛ 236line 235 didn't jump to line 236 because the condition on line 235 was never true
236 raise ExpressionParsingError(f"Cannot add constraint for {component}: only time-indexed components are supported.")
237 try:
238 def constraint_rule(blk, time_index):
239 expression = parse_expression(expr_str, fs,time_index)
240 # make sure the units of the expression are the same as the component
241 u1, u2 = get_attached_unit(component), get_attached_unit(expression)
242 if not check_units_equivalent(u1, u2): 242 ↛ 243line 242 didn't jump to line 243 because the condition on line 242 was never true
243 raise ValueError(
244 f"Failed to add constraint for {component}: units do not match (expected {u1}, got {u2})"
245 )
246 return pyunits.convert(component[time_index], to_units=u2) == expression
247 c = Constraint(component.index_set(), rule= constraint_rule)
248 name = f"equality_constraint_{id}"
249 fs.add_component(name, c)
250 add_corresponding_constraint(fs, c, id)
251 except ExpressionParsingError as e:
252 raise ExpressionParsingError(f"Failed to parse constraint expression '{expr_str}' for {component}: {expr_str}, error: {e}")
255 def build_spec_dependency_tree(self, specs) -> tuple:
256 """
257 Builds dependency tree for expressions based on the references in their respective expressions.
258 """
259 dependencies = defaultdict(
260 set
261 ) # Maps an expression's result_id to a set of result_ids it depends on
262 result_ids = set() # A set to track all result_ids
264 # Get a list of all result_id's.
265 for spec in specs:
266 result_id = spec["id"]
267 result_ids.add(result_id)
269 for spec in specs:
270 result_id = spec["id"]
271 expression = spec["expression"]
273 # Find the result_ids that this expression depends on
274 dependent_expressions = self.get_dependent_expressions(
275 expression, result_ids
276 )
278 if dependent_expressions:
279 # If the expression depends on another result_id, add dependency
280 for id in dependent_expressions:
281 dependencies[id].add(result_id)
283 return dependencies, result_ids
285 def get_dependent_expressions(self, expression: str, all_result_ids: set) -> list:
286 """
287 Gets all result_ids referenced in the expression.
288 """
289 # match result_ids starting with 'id_' followed by numbers
290 ids = re.findall(r"\b(id_\d+)\b", expression)
292 # Filter the referenced_ids to only include those that are in all_result_ids
293 valid_referenced_ids = []
294 for id in ids:
295 # get the numeric part after "id_" and check if it's in the all_result_ids
296 numeric_id = int(id[3:])
297 if numeric_id in all_result_ids:
298 valid_referenced_ids.append(numeric_id)
300 return valid_referenced_ids
302 def topological_sort(self, dependencies: dict, result_ids: set) -> list:
303 """
304 Performs topological sorting on the specification dependency tree.
305 """
306 # Track teh in-degree count for all expressions (edges coming into it)
307 in_degree = defaultdict(int)
309 # Count dependencies for each expression
310 for result_id in result_ids:
311 for dep in dependencies[result_id]:
312 in_degree[dep] += 1
314 # Initialise the queue with result_ids that have no dependencies (in-degree 0)
315 dequeue = deque(
316 [result_id for result_id in result_ids if in_degree[result_id] == 0]
317 )
319 sorted_result_ids = []
321 while dequeue:
322 result_id = dequeue.popleft()
323 sorted_result_ids.append(result_id)
325 # loop thtough and decrement each dependent expression's in_degree, and append it to the deque if it has an in-degree of 0.
326 for dependent_result_id in dependencies[result_id]:
327 in_degree[dependent_result_id] -= 1
328 if in_degree[dependent_result_id] == 0:
329 dequeue.append(dependent_result_id)
331 # If there are any result_ids left with non-zero in-degree, a cycle exists so error.
332 if len(sorted_result_ids) != len(result_ids): 332 ↛ 333line 332 didn't jump to line 333 because the condition on line 332 was never true
333 raise ValueError(
334 "Cycle detected in the dependency graph. Check an expression does not reference itself!"
335 )
337 return sorted_result_ids
339 def initialise(self) -> None:
340 """
341 Expands the arcs and initialises the model
342 """
343 # check if initialisation is disabled for this scenario
344 if getattr(self.schema, "disable_initialization", False):
345 # if disable initialisation is set to True, then we don't need to initialise the model
346 _log.info("Initialisation is disabled for this scenario.")
348 # We need to "expand_arcs" to make them a bidirection link that actually imposes constraints on the model.
349 TransformationFactory("network.expand_arcs").apply_to(self.model)
351 self.timing.step_into("initialise_model")
353 # load tear guesses (including var/constraint unfixing & deactivation and/or equality constraint deactivation)
354 self.tears.load()
356 tears = self.tears._tears
358 def init_unit(unit):
359 if getattr(self.schema, "disable_initialization", False): 359 ↛ 360line 359 didn't jump to line 360 because the condition on line 359 was never true
360 return
361 if (
362 getattr(
363 self.schema,
364 "skip_initialization_for_units_with_initial_values",
365 False,
366 )
367 and self.unit_models.has_initial_values(unit)
368 ):
369 _log.info(
370 f"Skipping initialization for unit {unit} because initial values exist."
371 )
372 return
374 _log.info("Initializing unit %s", unit)
375 self.timing.add_timing(f"init_{unit.name}")
376 try:
377 #unit.display()
378 unit.initialize(outlvl=idaeslog.INFO)
379 #unit.report()
380 except InitializationError as e:
381 # A unit initializer is a warm-start procedure. Even when its
382 # temporary subproblem is locally infeasible, the final iterate
383 # can still improve the starting point for the full flowsheet.
384 _log.warning(
385 "Unit initialization did not converge for %s; continuing "
386 "with its current values: %s",
387 unit.name,
388 e,
389 exc_info=True,
390 )
391 except Exception as e:
392 details = (
393 f"Unit initialization failed for {unit.name}: "
394 f"{type(e).__name__}: {e}"
395 )
396 _log.exception(details)
397 raise RuntimeError(details) from e
399 self.timing.add_timing("setup_sequential_decomposition")
400 # Use SequentialDecomposition to initialise the model
401 seq = SequentialDecomposition(
402 run_first_pass=True,
403 )
404 seq.set_tear_set(tears)
405 # use create_graph to get the order of sequential decomposition, and also to
406 # find any units that are not connected to the sequential decomposition
407 G = seq.create_graph(self.model)
408 order = seq.calculation_order(G)
409 seq_blocks = []
410 for o in order:
411 # o is a list of "levels" of the graph,
412 # all items in the first level have to no dependencies, the next level only depends on the previous levels.
413 for b in o:
414 seq_blocks.append(b)
415 _log.info("Order of initialisation: %s", [blk.name for blk in seq_blocks])
416 # set all the tear guesses before running the decomposition
417 for arc in tears:
418 port = arc.destination
419 # guesses used are initial values for each var
420 guesses = {}
421 guesses = {key: get_value(var) for key, var in port.vars.items()}
422 _log.info("Guess for %s: %s", port, guesses)
424 self.timing.step_into("run_sequential_decomposition")
426 # sequential decomposition completes when all vars across port
427 # equalities are within tol of each other
428 seq.options["tol"] = 1e-2
429 seq.options["solve_tears"] = False
430 seq.run(self.model, init_unit)
432 self.timing.step_out()
433 self.timing.add_timing("initialise_disconnected_units")
434 # Initialise any unit model that is not connected to the sequential decomposition
435 for blk in self.model.fs.component_data_objects(
436 Block, descend_into=False, active=True
437 ):
438 ports = list(blk.component_objects(Port, descend_into=False))
439 if len(ports) == 0:
440 continue # if the block has no ports, then it is not a unit model
441 if blk in seq_blocks:
442 continue # already initialised by sequential decomposition
443 init_unit(blk)
445 # unfix guess vars
446 deactivate_fixed_guesses(self.model.fs.guess_vars)
448 self.timing.step_out()
450 def serialise(self) -> SolvedFlowsheetSchema:
451 self.timing.add_timing("serialise_model")
453 initial_values = {}
454 for unit_model_id, unit_model in self.unit_models._unit_models.items():
455 initial_values[str(unit_model_id)] = to_json(unit_model, return_dict=True, wts=StoreSpec.value())
457 solved_flowsheet = SolvedFlowsheetSchema(
458 group_id=self.schema.group_id,
459 properties=serialize_properties_map(self.model.fs),
460 initial_values=initial_values
461 )
463 return solved_flowsheet
465 def report_statistics(self) -> None:
466 """
467 Reports statistics about the model
468 """
469 report_statistics(self.model)
471 def diagnose_problems(self) -> None:
472 from idaes.core.util import DiagnosticsToolbox
474 _log.info("=== DIAGNOSTICS ===")
476 report_statistics(self.model)
477 dt = DiagnosticsToolbox(self.model)
478 dt.report_structural_issues()
479 dt.display_overconstrained_set()
480 dt.display_underconstrained_set()
481 #dt.display_components_with_inconsistent_units()
482 try:
483 dt.compute_infeasibility_explanation()
484 except Exception as e:
485 _log.warning(
486 "Error computing infeasibility explanation: %s",
487 e,
488 )
489 try:
490 dt.report_numerical_issues()
491 except Exception as e:
492 _log.warning("Error reporting numerical issues: %s", e)
493 try:
494 dt.display_near_parallel_constraints()
495 except Exception as e:
496 _log.warning("Error displaying near parallel constraints: %s", e)
498 dt.display_variables_at_or_outside_bounds()
499 dt.display_variables_near_bounds()
500 try:
501 dt.display_constraints_with_extreme_jacobians()
502 except Exception as e:
503 _log.warning("Error displaying constraints with extreme jacobians: %s", e)
504 dt.display_constraints_with_large_residuals()
505 _log.info("=== END DIAGNOSTICS ===")
507 def get_unit_diagnostics(self):
508 diagnostics: list[tuple[Component, str]] = []
509 for unit_model in self.unit_models._unit_models.values():
510 if hasattr(unit_model, "diagnose"):
511 try:
512 results : list[tuple[Component, str]] = unit_model.diagnose()
513 diagnostics.extend(results)
514 # The unit diagnostics methods may be unreliable, as each unit model has different methods.
515 # So we wrap in a try-except to avoid the whole diagnostics process failing if one unit model's diagnostics fails.
516 except Exception as e:
517 _log.warning("Error Diagnosing %s: %s", unit_model, e)
518 return diagnostics
520 def get_unit_diagnostics_report(self) -> list[UnitDiagnosticsFinding]:
521 """Gets the diagnostics, but finds the property value ID for each diagnostic based on the component it references."""
522 diagnostics = self.get_unit_diagnostics()
523 diagnostics_report = [
524 UnitDiagnosticsFinding(
525 severity="error",
526 property_value_id=self.properties_map.get_id_by_component(component),
527 message=message)
528 for component, message in diagnostics
529 ]
530 diagnostics_report.extend(
531 _get_watertap_solver_diagnostics(
532 solver_option=self.schema.solver_option,
533 unit_models=self.unit_models._unit_models.values(),
534 properties_map=self.properties_map,
535 )
536 )
537 diagnostics_report.extend(
538 [UnitDiagnosticsFinding(property_value_id=id,
539 message=f"High Infeasibility: {infeasibility:.4f}")
540 for id, infeasibility in
541 get_top_infeasibilities(self.properties_map, threshold=1e-3)
542 ]
543 )
545 return diagnostics_report
547 def degrees_of_freedom(self) -> int:
548 """
549 Returns the degrees of freedom of the model
550 """
551 return int(degrees_of_freedom(self.model))
553 def check_model_valid(self) -> None:
554 """
555 Checks if the model is valid by checking the
556 degrees of freedom. Will raise an exception if
557 the model is not valid.
558 """
559 self.timing.add_timing("check_model_valid")
561 degrees_of_freedom = self.degrees_of_freedom()
562 if degrees_of_freedom != 0: 562 ↛ 564line 562 didn't jump to line 564 because the condition on line 562 was never true
563 #self.model.display() # prints the vars/constraints for debugging
564 raise Exception(
565 f"Degrees of freedom is not 0. Degrees of freedom: {degrees_of_freedom}"
566 )
568 def solve(self):
569 """
570 Solves the model. Uses petsc_dae_by_time_element for the petsc solver
571 (dynamic only), or a standard SolverFactory solver otherwise.
572 """
573 self.timing.add_timing("solve_model")
574 _log.info("=== Starting Solve ===")
576 if self.schema.solver_option == "petsc": 576 ↛ 577line 576 didn't jump to line 577 because the condition on line 576 was never true
577 from idaes.core.solvers import petsc as idaes_petsc
579 if not self.schema.dynamic:
580 raise ValueError("The PETSc solver is only supported for dynamic flowsheets.")
581 result = idaes_petsc.petsc_dae_by_time_element(
582 self.model,
583 time=self.model.fs.time,
584 ts_options={
585 "--ts_type": "beuler",
586 "--ts_dt": 1,
587 "--ts_monitor": "",
588 "--ts_save_trajectory": 1,
589 },
590 )
592 opt = SolverFactory(self.schema.solver_option)
593 # opt.options["max_iter"] = 5000
595 opt.options["max_iter"] = 1000
596 try:
597 starting_values = get_current_values(self.properties_map)
598 res = opt.solve(self.model, tee=True)
599 updated_values = get_current_values(self.properties_map)
600 changes = detect_changes(updated_values, starting_values)
601 print_changes(changes, self.properties_map)
603 success = check_optimal_termination(res)
604 return success
605 except ValueError as e:
606 if str(e).startswith("No variables appear"): 606 ↛ 610line 606 didn't jump to line 610 because the condition on line 606 was always true
607 # https://github.com/Pyomo/pyomo/pull/3445
608 return True
609 else:
610 raise e
612 def optimize(self) -> None:
613 if self.schema.optimizations is None or self.schema.optimizations == []:
614 return
616 self.timing.add_timing("optimize_model")
617 _log.info("=== Starting Optimization ===")
619 # ipopt doesn't support multiple objectives, so we need to create
620 # a single objective expression.
621 # this is done by summing all objectives in the model, adding
622 # or subtracting based on the sense (minimize or maximize)
623 objective_terms = []
624 for schema in self.schema.optimizations:
625 # get the expression component to optimize
626 # TODO: This is assuming optimisation is run on a steady-state simulation. We need to change how this works to handle dynamics.
627 # For now, just hardcoding time_index=0
628 objective_component = self.model.fs.properties_map.get_component(schema.objective)
629 # the objective component should be a scalar in non-dynamic models
630 # in a dynamic model it'll be indexed across all time steps
631 # We sum up all time steps so each time step is weighted equally.
632 objective = dimensionless_objective_expression(
633 sum_pyomo_expressions(objective_component.values())
634 )
635 # add or subtract the objective based on the sense
636 sense = schema.sense
637 if sense == "minimize":
638 objective_terms.append(objective)
639 else:
640 objective_terms.append(-objective)
644 # unfix relevant vars (add to degrees of freedom)
645 for dof_info in schema.unfixed_variables:
646 id = dof_info.id # Id of the propertyValue for this degree of freedom
647 var: PropertyComponent = self.model.fs.properties_map.get(id)
650 # TODO: may need to handle deactivating constraints,
651 # for expressions that are constrained (instead of state vars)
652 c = var.corresponding_constraint
653 if c is not None:
654 # TODO: better typing for constraints
655 if isinstance(c, ScalarConstraint) or isinstance(c, IndexedConstraint): 655 ↛ 656line 655 didn't jump to line 656 because the condition on line 655 was never true
656 c.deactivate()
657 else:
658 deactivate_components(c)
659 else:
660 for i in var.component.values():
661 if isinstance(i, ScalarVar): 661 ↛ 662line 661 didn't jump to line 662 because the condition on line 661 was never true
662 i.unfix()
663 # Because if not, it is ExpressionData, meaning it is already an expression and doesn't need to be unfixed. (we've already checked if there is a constraint for it above too.)
666 # TODO: set attributes for upper and lower bounds of property infos. i.e. use propertyinfo id.
667 # Var is either a Variable or Expression
668 # set the minimum or maximum bounds for this variable if they are enabled
669 #self.model.upper_bound_12 = Constraint(expr= var <= upper_bound_value )
671 upper_bound = dof_info.upper_bound
672 lower_bound = dof_info.lower_bound
674 c = var.component
676 if upper_bound is not None: 676 ↛ 682line 676 didn't jump to line 682 because the condition on line 676 was always true
677 def upper_bound_rule(model,index):
678 return c[index] <= upper_bound
679 upper_bound_constraint = Constraint(c.index_set(),rule=upper_bound_rule)
680 setattr(self.model,"upper_bound_" + str(id), upper_bound_constraint)
682 if lower_bound is not None: 682 ↛ 645line 682 didn't jump to line 645 because the condition on line 682 was always true
683 def lower_bound_rule(model,index):
684 return c[index] >= lower_bound
685 lower_bound_constraint = Constraint(c.index_set(),rule=lower_bound_rule)
686 setattr(self.model,"lower_bound_" + str(id), lower_bound_constraint)
689 # add the objective to the model
690 objective_expr = sum_pyomo_expressions(objective_terms)
691 self.model.objective = Objective(expr=objective_expr, sense=minimize)
693 # solve the model with the objective
694 opt = SolverFactory(self.schema.solver_option)
696 if self.schema.solver_option != "conopt": 696 ↛ 699line 696 didn't jump to line 699 because the condition on line 696 was always true
697 opt.options["max_iter"] = 1000
699 try:
700 res = opt.solve(self.model, tee=True)
701 assert_optimal_termination(res)
702 except ValueError as e:
703 if str(e).startswith("No variables appear"):
704 # https://github.com/Pyomo/pyomo/pull/3445
705 pass
706 else:
707 raise e