diff --git a/.gitignore b/.gitignore index 329184d..40a4c2b 100644 --- a/.gitignore +++ b/.gitignore @@ -472,3 +472,4 @@ $RECYCLE.BIN/ # End of https://www.toptal.com/developers/gitignore/api/python,jupyternotebooks,vim,visualstudiocode,pycharm,emacs,linux,macos,windows boulder/version.py +mixed_reactor_stream.yaml diff --git a/boulder/cantera_converter.py b/boulder/cantera_converter.py index 5e841aa..880fd03 100644 --- a/boulder/cantera_converter.py +++ b/boulder/cantera_converter.py @@ -34,6 +34,7 @@ class BoulderPlugins: reactor_builders: Dict[str, ReactorBuilder] = field(default_factory=dict) connection_builders: Dict[str, ConnectionBuilder] = field(default_factory=dict) post_build_hooks: List[PostBuildHook] = field(default_factory=list) + mechanism_path_resolver: Optional[Callable[[str], str]] = None # Global cache to ensure plugins are discovered only once @@ -112,7 +113,14 @@ def __init__( self.mechanism = mechanism or CANTERA_MECHANISM self.plugins = plugins or get_plugins() try: - self.gas = ct.Solution(self.mechanism) + # Use plugin mechanism path resolver if available + if self.plugins.mechanism_path_resolver: + resolved_mechanism = self.plugins.mechanism_path_resolver( + self.mechanism + ) + else: + resolved_mechanism = self.mechanism + self.gas = ct.Solution(resolved_mechanism) except Exception as e: raise ValueError(f"Failed to load mechanism '{self.mechanism}': {e}") self.reactors: Dict[str, ct.Reactor] = {} @@ -137,8 +145,28 @@ def create_reactor(self, reactor_config: Dict[str, Any]) -> ct.Reactor: reactor_type = reactor_config["type"] props = reactor_config["properties"] + # Determine mechanism for this node (override allowed via properties.mechanism) + mechanism_override = props.get("mechanism") + if mechanism_override: + try: + # Use plugin mechanism path resolver if available + if self.plugins.mechanism_path_resolver: + resolved_mechanism = self.plugins.mechanism_path_resolver( + mechanism_override + ) + else: + resolved_mechanism = mechanism_override + gas_obj = ct.Solution(resolved_mechanism) + except Exception as e: # fail-fast with clear message + raise ValueError( + f"Failed to load per-node mechanism '{mechanism_override}' for node " + f"'{reactor_config.get('id', '')}': {e}" + ) from e + else: + gas_obj = self.gas + # Set gas state - self.gas.TPX = ( + gas_obj.TPX = ( props.get("temperature", 300), props.get("pressure", 101325), self.parse_composition(props.get("composition", "N2:1")), @@ -148,19 +176,27 @@ def create_reactor(self, reactor_config: Dict[str, Any]) -> ct.Reactor: if reactor_type in self.plugins.reactor_builders: reactor = self.plugins.reactor_builders[reactor_type](self, reactor_config) elif reactor_type == "IdealGasReactor": - reactor = ct.IdealGasReactor(self.gas) + reactor = ct.IdealGasReactor(gas_obj) + elif reactor_type == "IdealGasMoleReactor": + reactor = ct.IdealGasMoleReactor(gas_obj) elif reactor_type == "IdealGasConstPressureReactor": - reactor = ct.IdealGasConstPressureReactor(self.gas) + reactor = ct.IdealGasConstPressureReactor(gas_obj) elif reactor_type == "IdealGasConstPressureMoleReactor": - reactor = ct.IdealGasConstPressureMoleReactor(self.gas) + reactor = ct.IdealGasConstPressureMoleReactor(gas_obj) elif reactor_type == "Reservoir": - reactor = ct.Reservoir(self.gas) + reactor = ct.Reservoir(gas_obj) else: raise ValueError(f"Unsupported reactor type: {reactor_type}") # Set the reactor name to match the config ID reactor.name = reactor_config["id"] + # Track per-node mechanism for downstream tools (e.g., sim2stone) + try: + reactor._boulder_mechanism = mechanism_override or self.mechanism # type: ignore[attr-defined] + except Exception: + pass + # Optional grouping: propagate group name to reactor for downstream tools props_group = reactor_config.get("properties", {}).get("group") if props_group is None: @@ -201,9 +237,11 @@ def create_connection(self, conn_config: Dict[str, Any]): # Default MassFlowController implementation mfc = ct.MassFlowController(source, target) mfc.mass_flow_rate = float(props.get("mass_flow_rate", 0.1)) + flow_device = mfc elif conn_type == "Valve": valve = ct.Valve(source, target) valve.valve_coeff = float(props.get("valve_coeff", 1.0)) + flow_device = valve elif conn_type == "Wall": # Handle walls as energy connections (e.g., torch power or losses) # After validation, electric_power_kW is converted to kilowatts if it had units @@ -371,11 +409,18 @@ def __init__( """ # Use provided mechanism or fall back to config default self.mechanism = mechanism or CANTERA_MECHANISM + self.plugins = plugins or get_plugins() try: - self.gas = ct.Solution(self.mechanism) + # Use plugin mechanism path resolver if available + if self.plugins.mechanism_path_resolver: + resolved_mechanism = self.plugins.mechanism_path_resolver( + self.mechanism + ) + else: + resolved_mechanism = self.mechanism + self.gas = ct.Solution(resolved_mechanism) except Exception as e: raise ValueError(f"Failed to load mechanism '{self.mechanism}': {e}") - self.plugins = plugins or get_plugins() self.reactors: Dict[str, ct.Reactor] = {} self.reactor_meta: Dict[str, Dict[str, Any]] = {} self.connections: Dict[str, ct.FlowDevice] = {} @@ -447,6 +492,51 @@ def build_network_and_code( ) except Exception: pass + elif typ == "IdealGasConstPressureReactor": + self.code_lines.append(f"{rid} = ct.IdealGasConstPressureReactor(gas)") + self.code_lines.append(f"{rid}.name = '{rid}'") + self.reactors[rid] = ct.IdealGasConstPressureReactor(self.gas) + self.reactors[rid].name = rid + try: + self.reactors[rid].group_name = str( + props.get("group", props.get("group_name", "")) + ) + self.code_lines.append( + f"{rid}.group_name = '{props.get('group', props.get('group_name', ''))}'" + ) + except Exception: + pass + elif typ == "IdealGasConstPressureMoleReactor": + self.code_lines.append( + f"{rid} = ct.IdealGasConstPressureMoleReactor(gas)" + ) + self.code_lines.append(f"{rid}.name = '{rid}'") + self.reactors[rid] = ct.IdealGasConstPressureMoleReactor(self.gas) + self.reactors[rid].name = rid + try: + self.reactors[rid].group_name = str( + props.get("group", props.get("group_name", "")) + ) + self.code_lines.append( + f"{rid}.group_name = '{props.get('group', props.get('group_name', ''))}'" + ) + except Exception: + pass + elif typ == "IdealGasMoleReactor": + # Available in Cantera 3.x + self.code_lines.append(f"{rid} = ct.IdealGasMoleReactor(gas)") + self.code_lines.append(f"{rid}.name = '{rid}'") + self.reactors[rid] = ct.IdealGasMoleReactor(self.gas) # type: ignore[attr-defined] + self.reactors[rid].name = rid + try: + self.reactors[rid].group_name = str( + props.get("group", props.get("group_name", "")) + ) + self.code_lines.append( + f"{rid}.group_name = '{props.get('group', props.get('group_name', ''))}'" + ) + except Exception: + pass elif typ == "Reservoir": self.code_lines.append(f"{rid} = ct.Reservoir(gas)") self.code_lines.append(f"{rid}.name = '{rid}'") diff --git a/boulder/config.py b/boulder/config.py index 62fe35b..b392a84 100644 --- a/boulder/config.py +++ b/boulder/config.py @@ -135,7 +135,7 @@ def normalize_config(config: Dict[str, Any]) -> Dict[str, Any]: for node in normalized["nodes"]: if "type" not in node: # Find the type key (anything that's not id, metadata, etc.) - standard_fields = {"id", "metadata"} + standard_fields = {"id", "metadata", "mechanism"} type_keys = [k for k in node.keys() if k not in standard_fields] if type_keys: @@ -148,6 +148,11 @@ def normalize_config(config: Dict[str, Any]) -> Dict[str, Any]: node["properties"] = ( properties if isinstance(properties, dict) else {} ) + # If a node-level mechanism is present, move it into properties for internal use + if "mechanism" in node: + props = node.setdefault("properties", {}) + if "mechanism" not in props: + props["mechanism"] = node["mechanism"] # Normalize connections if "connections" in normalized: @@ -295,10 +300,11 @@ def convert_to_stone_format(config: dict) -> dict: for node in config["nodes"]: # Build node with id first, then type as key containing properties node_type = node.get("type", "IdealGasReactor") - stone_node = { - "id": node["id"], - node_type: node.get("properties", {}), - } + props = dict(node.get("properties", {}) or {}) + mech_override = props.pop("mechanism", None) + stone_node = {"id": node["id"], node_type: props} + if mech_override is not None: + stone_node["mechanism"] = mech_override stone_config["nodes"].append(stone_node) # Convert connections (same structure in new format) diff --git a/boulder/sim2stone.py b/boulder/sim2stone.py new file mode 100644 index 0000000..99f2179 --- /dev/null +++ b/boulder/sim2stone.py @@ -0,0 +1,340 @@ +"""Conversion utilities: Cantera ReactorNet -> STONE YAML configuration. + +This module inspects a resolved Cantera network (``ct.ReactorNet``) and emits a +configuration compatible with Boulder/STONE. The conversion attempts to be +round-trippable: building a network from the emitted YAML should yield an +equivalent topology (same nodes and connections) when reconstructed by +``CanteraConverter``. +""" + +from __future__ import annotations + +import os +from typing import Any, Dict, List, Optional, Set, Tuple + +import cantera as ct # type: ignore + +from .config import CANTERA_MECHANISM, yaml_to_string_with_comments +from .ctutils import collect_all_reactors_and_reservoirs + + +def _mechanism_from_thermo(thermo: ct.ThermoPhase) -> Optional[str]: + """Try to read the mechanism file name from the ThermoPhase. + + Prefer attributes exposed by Cantera (e.g., `source`, `input_name`). Returns + a basename like 'gri30.yaml' when possible. + """ + for attr in ("source", "input_name"): + val = getattr(thermo, attr, None) + if isinstance(val, str) and val: + base = os.path.basename(val) + return base or val + return None + + +def _composition_to_string(thermo: ct.ThermoPhase, min_fraction: float = 1e-12) -> str: + """Serialize the current mole fractions to a compact composition string. + + Parameters + ---------- + thermo + A Cantera ThermoPhase (e.g., ``reactor.thermo``). + min_fraction + Species with mole fraction below this threshold are omitted to keep the + output readable. + """ + species: List[str] = list(thermo.species_names) + X: List[float] = list(thermo.X) + + # Filter near-zero entries for readability while preserving normalization. + items: List[Tuple[str, float]] = [ + (name, float(val)) + for name, val in zip(species, X) + if float(val) >= min_fraction + ] + if not items: + # Fallback when all are below threshold + return "N2:1" + + # Sort by descending fraction for stability + items.sort(key=lambda it: it[1], reverse=True) + return ",".join([f"{name}:{value:g}" for name, value in items]) + + +def _infer_node_type(reactor: ct.Reactor) -> str: + """Map a Cantera reactor instance to a Boulder component type string.""" + # Explicit checks for common types we support out-of-the-box + if isinstance(reactor, ct.Reservoir): + return "Reservoir" + if isinstance(reactor, ct.IdealGasReactor): + return "IdealGasReactor" + # Fallback to class name for supported custom/reactor variants + return type(reactor).__name__ + + +def _infer_connection_type(device: ct.FlowDevice) -> str: + """Map a Cantera flow device to a Boulder connection type string.""" + if isinstance(device, ct.MassFlowController): + return "MassFlowController" + if isinstance(device, ct.Valve): + return "Valve" + return type(device).__name__ + + +def _unique_flow_devices(all_reactors: Set[ct.Reactor]) -> Set[ct.FlowDevice]: + """Collect a unique set of flow devices from all reactors in the network.""" + devices: Set[ct.FlowDevice] = set() + for r in all_reactors: + # Both inlets and outlets are FlowDevice instances + for dev in list(getattr(r, "inlets", [])) + list(getattr(r, "outlets", [])): + devices.add(dev) + return devices + + +def _unique_walls(all_reactors: Set[ct.Reactor]) -> Set[ct.Wall]: + """Collect a unique set of walls from all reactors in the network.""" + walls: Set[ct.Wall] = set() + for r in all_reactors: + for w in getattr(r, "walls", []): + walls.add(w) + return walls + + +def sim_to_internal_config( + sim: ct.ReactorNet, default_mechanism: Optional[str] = None +) -> Dict[str, Any]: + """Convert a Cantera ReactorNet to Boulder internal configuration format. + + The returned config uses the normalized internal format consumed by + ``CanteraConverter.build_network``. + """ + mechanism = default_mechanism or CANTERA_MECHANISM + + # Collect unique reactors/reservoirs from the entire network + all_reactors = list(collect_all_reactors_and_reservoirs(sim)) + # Deterministic order + all_reactors.sort(key=lambda r: r.name) + + nodes: List[Dict[str, Any]] = [] + for r in all_reactors: + node_type = _infer_node_type(r) + # Capture mechanism override if present on reactor object (set by builder) + mech_override = getattr(r, "_boulder_mechanism", None) + + props: Dict[str, Any] = { + "temperature": float(r.thermo.T), + "pressure": float(r.thermo.P), + "composition": _composition_to_string(r.thermo), + } + if isinstance(mech_override, str) and mech_override: + props["mechanism"] = mech_override + else: + guessed = _mechanism_from_thermo(r.thermo) + if guessed: + props["mechanism"] = guessed + # Optional group propagation if present + group_name = getattr(r, "group_name", "") + if isinstance(group_name, str) and group_name: + props["group"] = group_name + + # Ensure each node has a non-empty unique identifier + rid = r.name or f"reactor_{id(r)}" + nodes.append({"id": rid, "type": node_type, "properties": props}) + + # Flow devices (MassFlowController, Valve, etc.) + devices = _unique_flow_devices(set(all_reactors)) + connections: List[Dict[str, Any]] = [] + for dev in sorted(list(devices), key=lambda d: (_infer_connection_type(d), id(d))): + src = dev.upstream + tgt = dev.downstream + src_id = getattr(src, "name", None) or f"reactor_{id(src)}" + tgt_id = getattr(tgt, "name", None) or f"reactor_{id(tgt)}" + + conn_type = _infer_connection_type(dev) + conn_props: Dict[str, Any] = {} + if isinstance(dev, ct.MassFlowController): + # Prefer attribute if available; Cantera exposes mass_flow_rate as a property + try: + props["mass_flow_rate"] = float(dev.mass_flow_rate) + except Exception: + # Fallbacks: some backends require initialized networks; try alternate attributes + mdot_attr = getattr(dev, "mdot", None) + try: + mdot_value = mdot_attr() if callable(mdot_attr) else mdot_attr + except Exception: + mdot_value = None + if isinstance(mdot_value, (int, float)): + props["mass_flow_rate"] = float(mdot_value) + # Else: omit property; builder will use its default + # If mass flow is negative, re-orient the connection and take absolute value + mfr = conn_props.get("mass_flow_rate") + if isinstance(mfr, (int, float)) and mfr < 0: + conn_props["mass_flow_rate"] = abs(float(mfr)) + # swap + src_id, tgt_id = tgt_id, src_id + elif isinstance(dev, ct.Valve): + # Capture valve coefficient if available + coeff = None + if hasattr(dev, "valve_coeff"): + coeff = getattr(dev, "valve_coeff") + elif hasattr(dev, "K"): + coeff = getattr(dev, "K") + if coeff is not None: + try: + conn_props["valve_coeff"] = float(coeff) + except Exception: + pass + + # Create a stable identifier + name_attr = getattr(dev, "name", None) + cid = ( + name_attr + if isinstance(name_attr, str) and name_attr + else f"{conn_type}_{src_id}_to_{tgt_id}" + ) + + connections.append( + { + "id": cid, + "type": conn_type, + "properties": conn_props, + "source": src_id, + "target": tgt_id, + } + ) + + # Walls (energy links) + walls = _unique_walls(set(all_reactors)) + for w in sorted( + list(walls), key=lambda w: (id(w.left_reactor), id(w.right_reactor)) + ): + left = w.left_reactor + right = w.right_reactor + l_id = getattr(left, "name", None) or f"reactor_{id(left)}" + r_id = getattr(right, "name", None) or f"reactor_{id(right)}" + + cid_raw = getattr(w, "name", None) + if not isinstance(cid_raw, str) or not cid_raw: + cid = f"Wall_{l_id}_to_{r_id}" + else: + cid = cid_raw + + # Convert current heat rate to an equivalent electric power in kW (best effort) + try: + q_watts = float(getattr(w, "heat_rate")) + except Exception: + q_watts = 0.0 + # If heat flows from right to left (negative), invert orientation and use positive magnitude + if q_watts < 0: + electric_power_kW = (-q_watts) / 1e3 + l_id, r_id = r_id, l_id + else: + electric_power_kW = q_watts / 1e3 + + connections.append( + { + "id": cid, + "type": "Wall", + "properties": { + "electric_power_kW": electric_power_kW, + # Efficiency unknown; preserve neutral defaults used by builder + "torch_eff": 1.0, + "gen_eff": 1.0, + }, + "source": l_id, + "target": r_id, + } + ) + + internal: Dict[str, Any] = { + "nodes": nodes, + "connections": connections, + "simulation": {"mechanism": mechanism}, + } + + return internal + + +def sim_to_stone_yaml( + sim: ct.ReactorNet, + default_mechanism: Optional[str] = None, +) -> str: + """Convert a Cantera ReactorNet to a STONE YAML string. + + Includes top-level `phases` (with `gas/mechanism`) and an explicit `settings` section. + Adds a blank line before `nodes` and `connections` for readability. + """ + internal = sim_to_internal_config(sim, default_mechanism=default_mechanism) + + # Determine mechanism for phases/gas + mechanism = ( + (internal.get("simulation") or {}).get("mechanism") + or default_mechanism + or CANTERA_MECHANISM + ) + + # Build STONE format using ruamel structures to control formatting + from ruamel.yaml.comments import CommentedMap, CommentedSeq + + stone_cm = CommentedMap() + # phases + phases_cm = CommentedMap() + gas_cm = CommentedMap() + gas_cm["mechanism"] = mechanism + phases_cm["gas"] = gas_cm + stone_cm["phases"] = phases_cm + + # settings (explicit section even if empty) + stone_cm["settings"] = CommentedMap() + + # nodes + nodes_seq = CommentedSeq() + for node in internal.get("nodes", []): + node_cm = CommentedMap() + node_cm["id"] = node["id"] + # Copy properties to avoid mutating internal + props = dict(node.get("properties", {}) or {}) + # Extract per-node mechanism and emit at node level (not inside class block) + node_mech = props.pop("mechanism", None) + # Drop per-node mechanism if it matches the global phases gas mechanism + if node_mech == mechanism: + node_mech = None + node_cm[node["type"]] = props + if node_mech is not None: + node_cm["mechanism"] = node_mech + nodes_seq.append(node_cm) + stone_cm["nodes"] = nodes_seq + # blank line before nodes + try: + stone_cm.yaml_set_comment_before_after_key("nodes", before="\n") + except Exception: + pass + + # connections + conns_seq = CommentedSeq() + for conn in internal.get("connections", []): + conn_cm = CommentedMap() + conn_cm["id"] = conn["id"] + conn_cm[conn["type"]] = conn.get("properties", {}) + conn_cm["source"] = conn["source"] + conn_cm["target"] = conn["target"] + conns_seq.append(conn_cm) + stone_cm["connections"] = conns_seq + # blank line before connections + try: + stone_cm.yaml_set_comment_before_after_key("connections", before="\n") + except Exception: + pass + + return yaml_to_string_with_comments(stone_cm) + + +def write_sim_as_yaml( + sim: ct.ReactorNet, + file_path: str, + default_mechanism: Optional[str] = None, +) -> None: + """Serialize a Cantera ReactorNet to a STONE YAML file.""" + yaml_str = sim_to_stone_yaml(sim, default_mechanism=default_mechanism) + with open(file_path, "w", encoding="utf-8") as f: + f.write(yaml_str) diff --git a/boulder/sim2stone_cli.py b/boulder/sim2stone_cli.py new file mode 100644 index 0000000..4c280e8 --- /dev/null +++ b/boulder/sim2stone_cli.py @@ -0,0 +1,160 @@ +"""CLI to convert a Cantera simulation (Python script) to STONE YAML. + +This command executes a Python file, searches its global variables for exactly +one ``ct.ReactorNet`` instance, and emits a STONE YAML file describing the +network topology and initial states. +""" + +from __future__ import annotations + +import argparse +import os +import sys +from typing import Optional + +import cantera as ct # type: ignore + +from .sim2stone import write_sim_as_yaml + + +def parse_args(argv: list[str] | None = None) -> argparse.Namespace: + parser = argparse.ArgumentParser( + prog="sim2stone", + description=( + "Execute a Python file, find a single Cantera ReactorNet, and convert " + "it to a STONE YAML file." + ), + ) + parser.add_argument( + "input", + help=( + "Path to a Python file that builds a Cantera ct.ReactorNet. The script " + "is executed in an isolated namespace. One and only one ReactorNet must " + "be present among its globals." + ), + ) + parser.add_argument( + "-o", + "--output", + default=None, + help=( + "Path to output YAML file. Defaults to replacing .py with .yaml next to " + "the input file." + ), + ) + parser.add_argument( + "--mechanism", + default=None, + help=( + "Override mechanism to record in phases/gas. Defaults to internal " + "configuration or library default." + ), + ) + parser.add_argument( + "--var", + default=None, + help=( + "Name of the global variable holding the ct.ReactorNet in the input script. " + "If provided, selection is forced even if multiple networks exist." + ), + ) + parser.add_argument( + "--verbose", action="store_true", help="Enable verbose progress messages" + ) + return parser.parse_args(argv) + + +def _execute_and_find_network( + script_path: str, var_name: Optional[str] = None +) -> ct.ReactorNet: + """Execute a Python script and return the single ReactorNet it defines. + + Adds the script directory to sys.path to resolve relative imports. + """ + import runpy + + script_abspath = os.path.abspath(script_path) + if not os.path.isfile(script_abspath): + raise FileNotFoundError(f"Input file does not exist: {script_path}") + + script_dir = os.path.dirname(script_abspath) + if script_dir and script_dir not in sys.path: + sys.path.insert(0, script_dir) + + # Execute script in its own globals namespace + globals_dict = runpy.run_path(script_abspath, run_name="__main__") + + # If a specific variable name is provided, use it + if var_name: + if var_name not in globals_dict or not isinstance( + globals_dict[var_name], ct.ReactorNet + ): + available = [ + k for k, v in globals_dict.items() if isinstance(v, ct.ReactorNet) + ] + raise RuntimeError( + f"Variable '{var_name}' is not a ct.ReactorNet in script globals. " + f"Available ReactorNet variables: {', '.join(available) if available else 'none'}" + ) + return globals_dict[var_name] + + # Collect ReactorNet instances from globals + networks = [v for v in globals_dict.values() if isinstance(v, ct.ReactorNet)] + + if len(networks) == 0: + raise RuntimeError( + "No ct.ReactorNet object found in script globals. Please ensure the " + "script assigns the network to a global variable (e.g., 'sim')." + ) + if len(networks) > 1: + names = [k for k, v in globals_dict.items() if isinstance(v, ct.ReactorNet)] + raise RuntimeError( + "Multiple ct.ReactorNet objects found in script globals: " + + ", ".join(names) + ) + return networks[0] + + +def main(argv: Optional[list[str]] = None) -> int: + args = parse_args(argv) + + if args.verbose: + print(f"[sim2stone] Executing: {args.input}") + + network = _execute_and_find_network(args.input, var_name=args.var) + + # Determine output path + if args.output: + output_path = args.output + else: + base, _ = os.path.splitext(os.path.abspath(args.input)) + output_path = base + ".yaml" + + if args.verbose: + print( + f"[sim2stone] Found ReactorNet with {len(network.reactors)} reactor(s). " + f"Writing: {output_path}" + ) + + write_sim_as_yaml(network, output_path, default_mechanism=args.mechanism) + + # Always print the created file path so callers can capture/chain it + print(f"🪨 STONE file created in {output_path}") + + # Validate the generated YAML (fail fast on errors) + from .config import load_config_file, normalize_config, validate_config + + cfg = load_config_file(output_path) + normalized = normalize_config(cfg) + validated = validate_config(normalized) + num_nodes = len(validated.get("nodes", [])) + num_conns = len(validated.get("connections", [])) + print(f"Validated STONE YAML ({num_nodes} nodes, {num_conns} connections)") + + if args.verbose: + print(f"[sim2stone] Done: {output_path}") + return 0 + + +if __name__ == "__main__": + raise SystemExit(main(sys.argv[1:])) diff --git a/boulder/styles.py b/boulder/styles.py index 2501ba4..353c6a4 100644 --- a/boulder/styles.py +++ b/boulder/styles.py @@ -60,10 +60,28 @@ "target-arrow-shape": "triangle", "target-arrow-color": "#555", "line-color": "#555", + # Ensure we can control draw order of edges + "z-index-compare": "manual", + # Keep default edges below walls + "z-index": 5, + # Slightly wider so they remain visible under walls + "width": 6, "text-wrap": "wrap", "text-max-width": "80px", }, }, + { + # Emphasize Walls with a distinct color visible in light mode + "selector": "edge[type = 'Wall']", + "style": { + "line-color": "#D0021B", + "target-arrow-color": "#D0021B", + # Draw walls on top + "z-index": 20, + # Slightly narrower so default edges peek around them + "width": 4, + }, + }, ] # Dark theme cytoscape stylesheet @@ -123,9 +141,27 @@ "target-arrow-shape": "triangle", "target-arrow-color": "#ccc", "line-color": "#ccc", + # Ensure we can control draw order of edges + "z-index-compare": "manual", + # Keep default edges below walls + "z-index": 5, "text-wrap": "wrap", "text-max-width": "80px", "color": "#fff", + # Slightly wider so they remain visible under walls + "width": 6, + }, + }, + { + # Emphasize Walls with a distinct color visible in dark mode + "selector": "edge[type = 'Wall']", + "style": { + "line-color": "#FF4D4D", + "target-arrow-color": "#FF4D4D", + # Draw walls on top + "z-index": 20, + # Slightly narrower so default edges peek around them + "width": 3, }, }, ] diff --git a/examples/mixed_reactor_stream_example.py b/examples/mixed_reactor_stream_example.py new file mode 100644 index 0000000..41671ff --- /dev/null +++ b/examples/mixed_reactor_stream_example.py @@ -0,0 +1,133 @@ +"""Mixing two streams example. + +The Cantera example https://cantera.org/3.1/examples/python/reactors/mix1.html, +used here to be converted to a 🪨 STONE serialized .yaml format with +:py:func:`boulder.sim2stone.sim_to_stone_yaml` + +Since reactors can have multiple inlets and outlets, they can be used to +implement mixers, splitters, etc. In this example, air and methane are mixed +in stoichiometric proportions. Due to the low temperature, no reactions occur. +Note that the air stream and the methane stream use different reaction +mechanisms, with different numbers of species and reactions. When gas flows +from one reactor or reservoir to another one with a different reaction +mechanism, species are matched by name. If the upstream reactor contains a +species that is not present in the downstream reaction mechanism, it will be +ignored. In general, reaction mechanisms for downstream reactors should +contain all species that might be present in any upstream reactor. +""" + +import cantera as ct + +# For regeneration (see end of file) +from boulder.cantera_converter import CanteraConverter +from boulder.config import load_config_file, normalize_config, validate_config + +# %% +# Boulder addition +# ---------------- +# For conversion: +from boulder.sim2stone import write_sim_as_yaml + +# %% +# Set up the reactor network +# -------------------------- +# +# Use air for stream a. +gas_a = ct.Solution("air.yaml") +gas_a.TPX = 300.0, ct.one_atm, "O2:0.21, N2:0.78, AR:0.01" +rho_a = gas_a.density + +# %% +# Use GRI-Mech 3.0 for stream b (methane) and for the mixer. If it is desired +# to have a pure mixer, with no chemistry, use instead a reaction mechanism +# for gas_b that has no reactions. +gas_b = ct.Solution("gri30.yaml") +gas_b.TPX = 300.0, ct.one_atm, "CH4:1" +rho_b = gas_b.density + +# %% +# Create reservoirs for the two inlet streams and for the outlet stream. The +# upstream reservoirs could be replaced by reactors, which might themselves be +# connected to reactors further upstream. The outlet reservoir could be +# replaced with a reactor with no outlet, if it is desired to integrate the +# composition leaving the mixer in time, or by an arbitrary network of +# downstream reactors. +res_a = ct.Reservoir(gas_a, name="Air Reservoir") +res_b = ct.Reservoir(gas_b, name="Fuel Reservoir") +downstream = ct.Reservoir(gas_a, name="Outlet Reservoir") + +# %% +# Create a reactor for the mixer. A reactor is required instead of a +# reservoir, since the state will change with time if the inlet mass flow +# rates change or if there is chemistry occurring. +gas_b.TPX = 300.0, ct.one_atm, "O2:0.21, N2:0.78, AR:0.01" +mixer = ct.IdealGasReactor(gas_b, name="Mixer") + +# %% +# Create two mass flow controllers connecting the upstream reservoirs to the +# mixer, and set their mass flow rates to values corresponding to +# stoichiometric combustion. +mfc1 = ct.MassFlowController(res_a, mixer, mdot=rho_a * 2.5 / 0.21, name="Air Inlet") +mfc2 = ct.MassFlowController(res_b, mixer, mdot=rho_b * 1.0, name="Fuel Inlet") + +# %% +# Connect the mixer to the downstream reservoir with a valve. +outlet = ct.Valve(mixer, downstream, K=10.0, name="Valve") + +sim = ct.ReactorNet([mixer]) + +# %% +# Get the mixed state +# ------------------- +# +# Since the mixer is a reactor, we need to integrate in time to reach steady state. +sim.advance_to_steady_state() + +# view the state of the gas in the mixer +print(mixer.thermo.report()) + +# %% +# Show the network structure +# -------------------------- +try: + diagram = sim.draw(print_state=True, species="X") +except ImportError as err: + print(f"Unable to show network structure:\n{err}") + +# %% +# Boulder additions +# ----------------- + +# Serialize the network to a STONE YAML file +output_yaml = "mixed_reactor_stream.yaml" +write_sim_as_yaml(sim, output_yaml, default_mechanism="gri30.yaml") +print(f"Wrote STONE YAML to {output_yaml}") + +# ✨ We now have a STONE YAML representation of the simulation! + +# %% +# The other way around: +# --------------------- +# Regenerate the simulation from the STONE YAML file + + +# Write the simulation as a STONE YAML file +write_sim_as_yaml(sim, "mixed_reactor_stream.yaml", default_mechanism="gri30.yaml") + +# Now regenerate the network from the YAML file, and verify it is runnable +loaded = load_config_file("mixed_reactor_stream.yaml") +normalized = normalize_config(loaded) +validated = validate_config(normalized) +converter = CanteraConverter( + mechanism=validated.get("simulation", {}).get("mechanism", "gri30.yaml") +) +network, results = converter.build_network(validated) +print( + f"Rebuilt network with {len(converter.reactors)} nodes and " + f"{len(converter.connections)} connections." +) + +# %% Compare the new network to the original +expected_nodes = {res_a.name, res_b.name, mixer.name, downstream.name} +assert set(converter.reactors.keys()) == expected_nodes +assert len(converter.connections) == 3 diff --git a/pyproject.toml b/pyproject.toml index dbf112d..d83b055 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -36,6 +36,7 @@ dev = [ [project.scripts] boulder = "boulder.cli:main" +sim2stone = "boulder.sim2stone_cli:main" [tool.coverage.run] branch = true diff --git a/tests/test_sim2stone_mechanisms.py b/tests/test_sim2stone_mechanisms.py new file mode 100644 index 0000000..76cc326 --- /dev/null +++ b/tests/test_sim2stone_mechanisms.py @@ -0,0 +1,51 @@ +from __future__ import annotations + +import cantera as ct # type: ignore + +from boulder.config import load_yaml_string_with_comments, normalize_config +from boulder.sim2stone import sim_to_stone_yaml + + +def test_sim2stone_preserves_per_node_mechanisms() -> None: + """Two reactors using different mechanisms should be preserved in YAML. + + R1 uses air.yaml; R2 uses gri30.yaml. The emitted STONE YAML should carry + node-level mechanism overrides that survive normalization. + """ + gas_air = ct.Solution("air.yaml") + gas_gri = ct.Solution("gri30.yaml") + + gas_air.TPX = 300.0, ct.one_atm, "N2:0.78,O2:0.21,AR:0.01" + gas_gri.TPX = 300.0, ct.one_atm, "CH4:1" + + r1 = ct.IdealGasReactor(gas_air, name="R1") + r2 = ct.IdealGasReactor(gas_gri, name="R2") + + # Hint serializer if supported + try: + r1._boulder_mechanism = "air.yaml" # type: ignore[attr-defined] + except Exception: + pass + try: + r2._boulder_mechanism = "gri30.yaml" # type: ignore[attr-defined] + except Exception: + pass + + sim = ct.ReactorNet([r1, r2]) + + # Use default equal to air.yaml so only R2 requires a per-node override + yaml_str = sim_to_stone_yaml(sim, default_mechanism="air.yaml") + loaded = load_yaml_string_with_comments(yaml_str) + normalized = normalize_config(loaded) + + mech_by_id = { + node["id"]: node["properties"].get("mechanism") for node in normalized["nodes"] + } + global_mech = normalized.get("simulation", {}).get("mechanism", None) + + # R1 matches global, so mechanism may be omitted at node level + assert global_mech == "air.yaml" + assert mech_by_id.get("R1") is None # because is equal to global "air.yaml" + + # R2 uses gri30, so mechanism override must be present at node level after normalization + assert mech_by_id.get("R2") == "gri30.yaml" diff --git a/tests/test_sim_to_yaml.py b/tests/test_sim_to_yaml.py new file mode 100644 index 0000000..66cc013 --- /dev/null +++ b/tests/test_sim_to_yaml.py @@ -0,0 +1,86 @@ +from __future__ import annotations + +import cantera as ct # type: ignore + +from boulder.cantera_converter import CanteraConverter +from boulder.config import ( + load_yaml_string_with_comments, + normalize_config, + validate_config, +) +from boulder.sim2stone import sim_to_internal_config, sim_to_stone_yaml + + +def _build_test_network(): + """Create a simple network: two reservoirs -> reactor -> reservoir.""" + gas = ct.Solution("gri30.yaml") + + # Upstream reservoirs with different compositions + gas.TPX = 300.0, ct.one_atm, "O2:0.21, N2:0.78, AR:0.01" + res_a = ct.Reservoir(gas, name="Air Reservoir") + + gas.TPX = 300.0, ct.one_atm, "CH4:1" + res_b = ct.Reservoir(gas, name="Fuel Reservoir") + + # Mixer reactor + gas.TPX = 300.0, ct.one_atm, "O2:0.21, N2:0.78, AR:0.01" + mixer = ct.IdealGasReactor(gas, name="Mixer") + + # Downstream sink + downstream = ct.Reservoir(gas, name="Outlet Reservoir") + + # Flow devices + mfc1 = ct.MassFlowController( + res_a, mixer, mdot=res_a.thermo.density * 2.5 / 0.21, name="Air Inlet" + ) + mfc2 = ct.MassFlowController( + res_b, mixer, mdot=res_b.thermo.density * 1.0, name="Fuel Inlet" + ) + valve = ct.Valve(mixer, downstream, K=10.0, name="Valve") + + # Only reactors go into ReactorNet (exclude pure reservoirs) + sim = ct.ReactorNet([mixer]) + return sim, {"mfc1": mfc1, "mfc2": mfc2, "valve": valve} + + +def test_roundtrip_sim_to_yaml_and_back(): + sim, _ = _build_test_network() + + # Convert sim -> internal config + internal = sim_to_internal_config(sim, default_mechanism="gri30.yaml") + + # Basic shape checks + assert ( + len(internal["nodes"]) == 4 + ) # Air Reservoir, Fuel Reservoir, Mixer, Outlet Reservoir + assert len(internal["connections"]) == 3 # two MFCs + one Valve + + # Serialize to STONE YAML string + yaml_str = sim_to_stone_yaml(sim, default_mechanism="gri30.yaml") + loaded_with_comments = load_yaml_string_with_comments(yaml_str) + normalized = normalize_config(loaded_with_comments) + validated = validate_config(normalized) + + # Rebuild via CanteraConverter + converter = CanteraConverter( + mechanism=validated.get("simulation", {}).get("mechanism", "gri30.yaml") + ) + network, results = converter.build_network(validated) + + # Node parity: same set of reactor IDs + original_ids = {n["id"] for n in internal["nodes"]} + rebuilt_ids = set(converter.reactors.keys()) + assert original_ids == rebuilt_ids + + # Connection parity (Flow devices only; Walls would be handled separately) + original_flow_ids = { + c["id"] + for c in internal["connections"] + if c["type"] in ("MassFlowController", "Valve") + } + rebuilt_flow_ids = set(converter.connections.keys()) + assert original_flow_ids == rebuilt_flow_ids + + # ReactorNet parity: one non-reservoir reactor called "Mixer" + assert len(network.reactors) == 1 + assert network.reactors[0].name == "Mixer"