Coverage for backend/ahuora-builder/src/ahuora_builder/methods/expression_parsing.py: 91%

1"""Parse factory-emitted expression strings into Pyomo expressions. 

2 

3Supported expression-language functions are intentionally narrow. The Django 

4factory emits lowercase ``min(arg, arg, ...)`` and ``max(arg, arg, ...)`` for 

5property-key aggregates, plus ``convert(expr, unit)`` when generated formulas 

6need an explicit unit basis. Uppercase or mixed-case variants remain invalid so 

7that user-facing aggregate names are resolved before expressions reach the 

8builder. 

9""" 

10 

11import re 

12from typing import Any 

13from pyomo.environ import Expr_if, Expression, Param 

14from pyomo.core.base.expression import ExpressionData 

15from pyomo.core.base.units_container import units as pyomo_units, _PyomoUnit 

16from idaes.core.util.math import smooth_min 

17from sympy import Symbol 

18from sympy.parsing.sympy_parser import parse_expr 

19from pyomo.core.expr.sympy_tools import sympy2pyomo_expression, PyomoSympyBimap 

20from pyomo.core.base.indexed_component import IndexedComponent 

21from ahuora_builder.properties_manager import PropertiesManager 

22from ahuora_builder_types.units import register_project_units 

23from .slice_manipulation import is_scalar_reference 

24class ExpressionParsingError(Exception): 

25 """Custom exception for errors during expression parsing.""" 

26 pass 

27 

28# add Ahuora-specific units to Pyomo's Pint registry 

29ureg = pyomo_units._pint_registry 

30register_project_units(ureg.define) 

31 

32# Builder function grammar: lowercase special function identifiers followed by 

33# a parenthesized comma-separated argument list. Uppercase or mixed-case calls 

34# are rejected explicitly instead of being handed to SymPy. 

35SPECIAL_FUNCTION_NAMES = {"min", "max", "convert"} 

36SPECIAL_FUNCTION_PATTERN = re.compile(r"\b(min|max|convert)\s*\(") 

37SPECIAL_FUNCTION_CASE_VARIANT_PATTERN = re.compile(r"\b([A-Za-z_][A-Za-z0-9_]*)\s*\(") 

38SMOOTH_MIN_EPSILON = 1e-6 

39 

40 

41def handle_special_chars(expr: str) -> str: 

42 # replace special characters so they can be parsed 

43 expr = expr.replace("^", "**") 

44 expr = expr.replace("$", "dollar") 

45 

46 return expr 

47 

48 

49def get_property_from_id(fs, property_id, time_index): 

50 """Return the Pyomo property data object referenced by a property id.""" 

51 

52 properties_map : PropertiesManager = fs.properties_map 

53 pyomo_object: IndexedComponent = properties_map.get_component(property_id) 

54 

55 if pyomo_object is None: 55 ↛ 56line 55 didn't jump to line 56 because the condition on line 55 was never true

56 raise ValueError(f"Symbol with id {property_id} not found in model") 

57 # check if this is a time-indexed var, and if so get the value at the given time index 

58 if is_scalar_reference(pyomo_object): 

59 # reference with index None 

60 return pyomo_object[None] 

61 elif pyomo_object.index_set() == fs.time: 61 ↛ 64line 61 didn't jump to line 64 because the condition on line 61 was always true

62 return pyomo_object[time_index] 

63 else: 

64 raise NotImplementedError("Only 0D and 1D time-indexed properties are supported in expressions") 

65 

66 

67def get_property_expression_for_symbol(fs, property_id, time_index): 

68 """Return the value expression parsing should substitute for an ``id_`` symbol.""" 

69 

70 value = get_property_from_id(fs, property_id, time_index) 

71 if isinstance(value, ExpressionData): 

72 return value.expr 

73 return value 

74 

75 

76def evaluate_symbol(fs, symbol: str,time_index) -> Any: 

77 if symbol.lower() == "time" or symbol.lower() == "t": 77 ↛ 78line 77 didn't jump to line 78 because the condition on line 77 was never true

78 return float(time_index) 

79 if symbol.startswith("id_"): 

80 # get the property from flowsheet properties_map 

81 id = int(symbol[3:]) 

82 return get_property_expression_for_symbol(fs, id, time_index) 

83 else: 

84 # assume its a unit, eg. "m" or "kg" 

85 # get the unit from pint, pyomo's units library 

86 ureg = pyomo_units._pint_registry 

87 pint_unit = getattr(ureg, symbol) 

88 pyomo_unit = _PyomoUnit(pint_unit, ureg) 

89 # We want people to write expressions such as (10 * W + 5 * kW). Pyomo doesn't natively support this, 

90 # so we can always convert to base units. 

91 if symbol == "delta_degC" or symbol == "delta_degF": 

92 # special case, because degC is not a base unit 

93 return 1 * pyomo_unit 

94 #return _PyomoUnit(ureg.delta_degC) 

95 elif symbol == "degC" or symbol == "degF": 95 ↛ 98line 95 didn't jump to line 98 because the condition on line 95 was never true

96 # throw an error (we do not support this, as it is unclear what to do) 

97 # https://pyomo.readthedocs.io/en/6.8.1/explanation/modeling/units.html 

98 raise ValueError(f"Use relative temperature units (delta_degC, delta_degF) or absolute temperature units (K, degF). Cannot use {symbol} as addition and multiplication is inconsistent on non-absolute units") 

99 scale_factor, base_units = ureg.get_base_units(pint_unit, check_nonmult=True) # TODO: handle degC etc. 

100 base_pyomo_unit = _PyomoUnit(base_units, ureg) 

101 return pyomo_units.convert( 1 * pyomo_unit, to_units=base_pyomo_unit) 

102 

103 

104def split_function_args(arguments: str) -> list[str]: 

105 """Split function arguments without treating nested commas as separators.""" 

106 

107 args = [] 

108 start = 0 

109 depth = 0 

110 for index, char in enumerate(arguments): 

111 if char == "(": 

112 depth += 1 

113 elif char == ")": 

114 depth -= 1 

115 if depth < 0: 115 ↛ 116line 115 didn't jump to line 116 because the condition on line 115 was never true

116 raise ExpressionParsingError( 

117 f"Unexpected closing parenthesis in function arguments '{arguments}'" 

118 ) 

119 elif char == "," and depth == 0: 

120 args.append(arguments[start:index].strip()) 

121 start = index + 1 

122 args.append(arguments[start:].strip()) 

123 if depth != 0: 123 ↛ 124line 123 didn't jump to line 124 because the condition on line 123 was never true

124 raise ExpressionParsingError( 

125 f"Unbalanced parentheses in function arguments '{arguments}'" 

126 ) 

127 return args 

128 

129 

130def find_matching_paren(expression: str, open_index: int) -> int: 

131 """Find the closing parenthesis matching the opening parenthesis at open_index.""" 

132 

133 depth = 0 

134 for index in range(open_index, len(expression)): 

135 if expression[index] == "(": 

136 depth += 1 

137 elif expression[index] == ")": 

138 depth -= 1 

139 if depth == 0: 

140 return index 

141 raise ExpressionParsingError(f"Unclosed function call in expression '{expression}'") 

142 

143 

144def _smooth_min_epsilon(model, expression: Any) -> Any: 

145 """Return a small smoothing term compatible with IDAES ``smooth_min``. 

146 

147 IDAES accepts a numeric epsilon for unitless expressions, but unitful 

148 expressions need a Pyomo ``Param`` carrying the same units so Pyomo can 

149 prove ``sqrt((a - b)^2 + eps^2)`` is dimensionally consistent. 

150 """ 

151 

152 expression_units = pyomo_units.get_units(expression) 

153 if expression_units is None: 153 ↛ 154line 153 didn't jump to line 154 because the condition on line 153 was never true

154 return SMOOTH_MIN_EPSILON 

155 

156 unit_key = str(expression_units) 

157 epsilon_cache = getattr(model, "_ahuora_smooth_min_epsilon_by_unit", None) 

158 if epsilon_cache is None: 

159 epsilon_cache = {} 

160 setattr(model, "_ahuora_smooth_min_epsilon_by_unit", epsilon_cache) 

161 if unit_key in epsilon_cache: 

162 return epsilon_cache[unit_key] 

163 

164 component_index = len(epsilon_cache) 

165 while True: 

166 component_name = f"_ahuora_smooth_min_epsilon_{component_index}" 

167 if not hasattr(model, component_name): 167 ↛ 169line 167 didn't jump to line 169 because the condition on line 167 was always true

168 break 

169 component_index += 1 

170 

171 epsilon = Param(initialize=SMOOTH_MIN_EPSILON, units=expression_units) 

172 model.add_component(component_name, epsilon) 

173 epsilon_cache[unit_key] = epsilon 

174 return epsilon 

175 

176 

177def _pairwise_min_max(model, function_name: str, args: list[Any]) -> Any: 

178 """Build a left-folded expression-language min/max expression. 

179 

180 ``min`` lowers to IDAES ``smooth_min`` so objective functions remain 

181 differentiable. ``max`` keeps the existing exact conditional semantics 

182 until there is a product requirement to smooth it too. 

183 """ 

184 

185 expression = args[0] 

186 for arg in args[1:]: 

187 expression_units = pyomo_units.get_units(expression) 

188 arg_units = pyomo_units.get_units(arg) 

189 if expression_units is not None and arg_units is not None: 189 ↛ 191line 189 didn't jump to line 191 because the condition on line 189 was always true

190 arg = pyomo_units.convert(arg, to_units=expression_units) 

191 if function_name == "min": 

192 expression = smooth_min(expression, arg, eps=_smooth_min_epsilon(model, expression)) 

193 else: 

194 expression = Expr_if(IF_=expression >= arg, THEN_=expression, ELSE_=arg) 

195 return expression 

196 

197 

198def _convert_to_units(expression: Any, target_quantity: Any) -> Any: 

199 """Convert an expression to the units carried by a parsed unit quantity.""" 

200 

201 target_units = pyomo_units.get_units(target_quantity) 

202 if target_units is None: 202 ↛ 203line 202 didn't jump to line 203 because the condition on line 202 was never true

203 raise ExpressionParsingError("convert target must include units") 

204 return pyomo_units.convert(expression, to_units=target_units) 

205 

206 

207def replace_special_functions(expression: str, model, time_index) -> tuple[str, dict[str, Any]]: 

208 """Lower expression-language special calls before handing the rest to SymPy. 

209 

210 SymPy can parse the names ``min`` and ``max``, but the result is not a Pyomo 

211 expression. Django emits lowercase builder expression-language calls for 

212 aggregate ``MIN``/``MAX`` and generated unit coercions; this pass lowers 

213 those calls into Pyomo placeholders whose expressions are returned by the 

214 bimap unchanged. 

215 """ 

216 

217 special_symbols = {} 

218 

219 while True: 

220 case_variant = next( 

221 ( 

222 match 

223 for match in SPECIAL_FUNCTION_CASE_VARIANT_PATTERN.finditer(expression) 

224 if match.group(1).lower() in SPECIAL_FUNCTION_NAMES 

225 and match.group(1) != match.group(1).lower() 

226 ), 

227 None, 

228 ) 

229 if case_variant is not None: 

230 raise ExpressionParsingError( 

231 f"{case_variant.group(1)} must be written as lowercase " 

232 f"{case_variant.group(1).lower()}" 

233 ) 

234 

235 match = SPECIAL_FUNCTION_PATTERN.search(expression) 

236 if match is None: 

237 return expression, special_symbols 

238 

239 open_index = match.end() - 1 

240 close_index = find_matching_paren(expression, open_index) 

241 function_name = match.group(1) 

242 args = split_function_args(expression[open_index + 1:close_index]) 

243 if any(arg == "" for arg in args): 

244 raise ExpressionParsingError(f"{function_name} arguments cannot be empty") 

245 if function_name in {"min", "max"} and len(args) < 2: 

246 raise ExpressionParsingError( 

247 f"{function_name} expects at least 2 arguments, got {len(args)}" 

248 ) 

249 if function_name == "convert" and len(args) != 2: 

250 raise ExpressionParsingError( 

251 f"convert expects exactly 2 arguments, got {len(args)}" 

252 ) 

253 parsed_args = [parse_expression(arg, model, time_index) for arg in args] 

254 placeholder = f"__ahuora_special_{len(special_symbols)}" 

255 if function_name == "convert": 

256 special_symbols[placeholder] = _convert_to_units(parsed_args[0], parsed_args[1]) 

257 else: 

258 special_symbols[placeholder] = _pairwise_min_max(model, function_name, parsed_args) 

259 expression = expression[:match.start()] + placeholder + expression[close_index + 1:] 

260 

261 

262class PyomoSympyMap(PyomoSympyBimap): 

263 """Map expression symbols to Pyomo objects during SymPy-to-Pyomo conversion.""" 

264 

265 def __init__(self, model,time_index, special_symbols: dict[str, Any] | None = None): 

266 self.model = model 

267 self.time_index = time_index 

268 self.special_symbols = special_symbols or {} 

269 

270 def getPyomoSymbol(self, sympy_object: Symbol, default=None): 

271 if not isinstance(sympy_object, Symbol): 

272 return None # It's not in pyomo, e.g a number or something 

273 if sympy_object.name in self.special_symbols: 

274 return self.special_symbols[sympy_object.name] 

275 return evaluate_symbol(self.model, sympy_object.name, self.time_index) 

276 

277 def getSympySymbol(self, pyomo_object, default=None): 

278 raise NotImplementedError( 

279 "getSympySymbol not implemented, because it shouldn't be needed" 

280 ) 

281 # we don't care, it only needs to go one way 

282 

283 def sympyVars(self): 

284 raise NotImplementedError( 

285 "sympyVars not implemented, because it shouldn't be needed" 

286 ) 

287 

288 

289def parse_expression(expression, model,time_index) -> Expression: 

290 """Parse an Ahuora expression string into a Pyomo expression.""" 

291 

292 # use the bimap to get the correct pyomo object for each symbol 

293 try: 

294 expression = handle_special_chars(expression) 

295 expression, special_symbols = replace_special_functions(expression, model, time_index) 

296 bimap = PyomoSympyMap(model,time_index, special_symbols) 

297 sympy_expr = parse_expr(expression) 

298 pyomo_expr = sympy2pyomo_expression(sympy_expr, bimap) 

299 except Exception as e: 

300 raise ExpressionParsingError(f"{expression}: error: {e}") 

301 return pyomo_expr