From 391ffa412062efb1517b112ade0fcba4b3192621 Mon Sep 17 00:00:00 2001 From: Yoel Cortes-Pena Date: Tue, 30 Jan 2024 00:34:52 -0600 Subject: [PATCH] improvement to PO-sim --- biosteam/_system.py | 33 ++++++++++++--------- biosteam/units/phase_equilibrium.py | 4 +-- tests/test_phenomena_oriented_simulation.py | 3 +- 3 files changed, 23 insertions(+), 17 deletions(-) diff --git a/biosteam/_system.py b/biosteam/_system.py index 7a14ebe4..7d76bbf1 100644 --- a/biosteam/_system.py +++ b/biosteam/_system.py @@ -120,24 +120,29 @@ def solve(self): values.append(value) return objs, np.array(values) A, objs = dictionaries2array(self.A) - b = np.array(b).T - if A.ndim == 2: + if A.ndim == 3: + A_ = A + b_ = np.array(b).T + objs_ = objs + values = [] + for A, b in zip(A_, b_): + rows = A.any(axis=0) + cols = A.any(axis=1) + b = [j for (i, j) in zip(rows, b)] + A = A[rows][:, cols] + objs = [j for (i, j) in zip(cols, objs_) if i] + values.append(solve(A, b).T) + values = np.array(values).T + else: rows = A.any(axis=0) cols = A.any(axis=1) + b = np.array(b).T b = [j for (i, j) in zip(rows, b)] A = A[rows][:, cols] objs = [j for (i, j) in zip(cols, objs) if i] - else: - # TODO: A.ndim == 3 - pass - try: values = solve(A, b).T - except Exception as e: - # for i in self.A: - # print('--') - # for i, j in i.items(): - # print(i.ID, j) - raise e + if np.isnan(values).any(): + raise RuntimeError('nan value in variables') for obj, value in zip(objs, values): obj._update_decoupled_variable(variable, value) if variable in ('mol', 'mol-LLE'): @@ -685,7 +690,7 @@ class System: available_methods: Methods[str, tuple[Callable, bool, dict]] = Methods() #: Variable solution priority for phenomena oriented simulation. - variable_priority: list[str] = ['mol', ('mol-LLE', 'K-pseudo'), 'K', 'T', 'L', 'B'] + variable_priority: list[str] = ['mol', ('mol-LLE', 'K-pseudo'), 'K', 'B', 'mol', 'T', 'L'] @classmethod def register_method(cls, name, solver, conditional=False, **kwargs): @@ -2222,7 +2227,7 @@ def run(self): self.run_phenomena() except: warn('phenomena-oriented simulation failed; ' - 'attempting one sequential-modular loop', RuntimeWarning) + 'attempting one sequential-modular loop', RuntimeWarning) for i in self.unit_path: i.run() else: raise RuntimeError(f'unknown algorithm {algorithm!r}') diff --git a/biosteam/units/phase_equilibrium.py b/biosteam/units/phase_equilibrium.py index c9a5d26f..81082e19 100644 --- a/biosteam/units/phase_equilibrium.py +++ b/biosteam/units/phase_equilibrium.py @@ -504,7 +504,7 @@ def _get_arrays(self): def _set_arrays(self, IDs, **kwargs): IDs_last = self.IDs - if IDs_last and IDs_last != IDs: + if IDs_last and IDs_last != IDs and len(IDs_last) > len(IDs): size = len(IDs_last) index = [IDs_last.index(i) for i in IDs] for name, array in kwargs.items(): @@ -668,7 +668,7 @@ def _run_vle(self, P=None, update=True): y_mol = 0 K_new = y_mol / x_mol if B is None: - if V_total and not L_total: + if not L_total: self.B = inf else: self.B = V_total / L_total diff --git a/tests/test_phenomena_oriented_simulation.py b/tests/test_phenomena_oriented_simulation.py index b29301b0..1807046d 100644 --- a/tests/test_phenomena_oriented_simulation.py +++ b/tests/test_phenomena_oriented_simulation.py @@ -197,7 +197,7 @@ def fresh_solvent_flow_rate(): for i in range(1): sm_sys.simulate() t1 = time.toc() - assert t0 < 0.75 * t1, 'phenomena-oriented simulation speed should be faster than sequential-modular' + assert t0 < 0.2 * t1, 'phenomena-oriented simulation speed should be faster than sequential-modular' for s_sm, s_dp in zip(sm_sys.streams, dp_sys.streams): actual = s_sm.mol @@ -296,6 +296,7 @@ def adjust_fresh_solvent_flow_rate(): time.tic() sm = create_system('sequential modular') + sm.flatten() sm.set_tolerance(rmol=1e-6, mol=1e-6, subsystems=True, method='fixed-point') sm.simulate() t_sequential = time.toc()