From cd91740d16484e9eadf060ed166d125e78bef2f9 Mon Sep 17 00:00:00 2001 From: Yoel Cortes-Pena Date: Thu, 30 May 2024 19:14:12 -0500 Subject: [PATCH] more enhancements to PO-API --- benchmark/acetic_acid_system.py | 10 +- benchmark/butanol_purification_system.py | 10 +- benchmark/profile.py | 19 +- biosteam/_system.py | 28 ++- biosteam/units/compressor.py | 4 +- biosteam/units/distillation.py | 8 + biosteam/units/mixing.py | 6 +- biosteam/units/splitting.py | 4 +- biosteam/units/stage.py | 194 +++++++++++++++----- tests/test_equilibrium_stages.py | 6 +- tests/test_phenomena_oriented_simulation.py | 2 +- 11 files changed, 205 insertions(+), 86 deletions(-) diff --git a/benchmark/acetic_acid_system.py b/benchmark/acetic_acid_system.py index 28f52375..ed454100 100644 --- a/benchmark/acetic_acid_system.py +++ b/benchmark/acetic_acid_system.py @@ -6,7 +6,11 @@ """ import biosteam as bst import numpy as np -from .profile import register +try: + from .profile import register +except: + def register(*args, **kwargs): + return lambda f: f __all__ = ( 'create_acetic_acid_simple_system', @@ -77,7 +81,7 @@ def fresh_solvent_flow_rate(): @register( 'acetic_acid_complex', 'Glacial acetic acid\npurification', - 80, [0, 10, 20, 30, 40, 50, 60, 70, 80], 'AcOH\nsep.' + 60, [0, 10, 20, 30, 40, 50, 60], 'AcOH\nsep.' ) def create_acetic_acid_complex_system(alg): thermo = bst.Thermo(['Water', 'AceticAcid', 'EthylAcetate'], cache=True) @@ -169,6 +173,7 @@ def adjust_fresh_solvent_flow_rate(): reflux=None, boilup=3, use_cache=True, + maxiter=10, ) settler = bst.StageEquilibrium( 'settler', @@ -226,6 +231,7 @@ def adjust_fresh_solvent_flow_rate(): feed_stages=(1, 2), reflux=1, boilup=2, + maxiter=10, ) return sys diff --git a/benchmark/butanol_purification_system.py b/benchmark/butanol_purification_system.py index 620bebe2..f06f0f70 100644 --- a/benchmark/butanol_purification_system.py +++ b/benchmark/butanol_purification_system.py @@ -5,7 +5,11 @@ @author: cortespea """ import biosteam as bst -from .profile import register +try: + from .profile import register +except: + def register(*args, **kwargs): + return lambda f: f __all__ = ( 'create_system_butanol_purification', @@ -13,7 +17,7 @@ @register( 'butanol_purification', 'Butanol purification', - 5, [1, 2, 3, 4, 5], 'BtOH\nsep.' + 2, [0.5, 1, 1.5, 2], 'BtOH\nsep.' ) def create_system_butanol_purification(alg): bst.settings.set_thermo(['Water', 'Butanol'], cache=True) @@ -53,5 +57,5 @@ def create_system_butanol_purification(alg): LHK=('Water', 'Butanol'), ) - sys = bst.System.from_units(units=[water_distiller, settler, butanol_distiller]) + sys = bst.System.from_units(units=[water_distiller, settler, butanol_distiller], algorithm=alg) return sys diff --git a/benchmark/profile.py b/benchmark/profile.py index 05c55339..772c969e 100644 --- a/benchmark/profile.py +++ b/benchmark/profile.py @@ -150,6 +150,8 @@ def profile_phenomena_oriented(system, total_time): po = all_systems[system]('phenomena oriented') po.flatten() po._setup_units() + # print(po.algorithm) + # breakpoint() data = profile(po.run_phenomena, *streams_and_stages(po), total_time) bst.F.clear() return data @@ -191,7 +193,6 @@ def dT_error(stage): sum([i.H for i in stage.outs]) - sum([i.H for i in stage.ins]) ) / sum([i.C for i in stage.outs]) -# @high_precision def benchmark(f, streams, adiabatic_stages, stages, total_time): time = bst.TicToc() net_time = 0 @@ -361,29 +362,29 @@ def plot_benchmark(systems=None, N=100, load=True, save=True): # @high_precision def profile(f, streams, adiabatic_stages, stages, total_time): time = bst.TicToc() - KB_error = [] flow_error = [] energy_error = [] material_error = [] temperature_error = [] record = [] net_time = 0 - KBs = np.array([i.K * i.B for i in stages if hasattr(i, 'K')]) temperatures = np.array([i.T for i in streams]) flows = np.array([i.mol for i in streams]) while net_time < total_time: time.tic() f() net_time += time.toc() - new_KBs = np.array([i.K * i.B for i in stages if hasattr(i, 'K')]) new_temperatures = np.array([i.T for i in streams]) new_flows = np.array([i.mol for i in streams]) record.append(net_time) dF = np.abs(flows - new_flows).sum() + # sys = f.__self__ + # print(sys.algorithm) + # if sys.algorithm == 'Phenomena oriented': + # sys.show() + # breakpoint() + # print(dF) dT = np.abs(temperatures - new_temperatures).sum() - KB_error.append( - np.log10(np.abs(new_KBs - KBs).sum() + 1e-15) - ) flow_error.append( np.log10(dF + 1e-15) ) @@ -396,12 +397,10 @@ def profile(f, streams, adiabatic_stages, stages, total_time): material_error.append( np.log10(sum([abs(i.mass_balance_error()) for i in stages]) + 1e-15) ) - KBs = new_KBs flows = new_flows temperatures = new_temperatures return { 'Time': record, - 'Stripping factor': KB_error, 'Stream temperature': temperature_error, 'Component flow rate': flow_error, 'Energy balance': energy_error, @@ -436,7 +435,7 @@ def dct_mean_std(dcts: list[dict], keys: list[str]): return {i: (values[i].mean(), values[i].std()) for i in keys} def plot_profile( - systems=None, N=10, load=True, save=True + systems=None, N=1, load=True, save=True ): if systems is None: systems = list(all_systems) fs = 9 diff --git a/biosteam/_system.py b/biosteam/_system.py index 7b8b70fc..368a2349 100644 --- a/biosteam/_system.py +++ b/biosteam/_system.py @@ -91,6 +91,8 @@ def solve_material_flows(self): A = [] b = [] for node in nodes: + if not node._create_material_balance_equations: + raise NotImplementedError(f'{node!r} has no method `_create_material_balance_equations`') for coefficients, value in node._create_material_balance_equations(): coefficients = { streams[i.imol]: j for i, j in coefficients.items() @@ -139,13 +141,17 @@ def solve_energy_departures(self): A = [] b = [] for node in nodes: + if not node._create_energy_departure_equations: + raise NotImplementedError(f'{node!r} has no method `_create_energy_departure_equations`') for coefficients, value in node._create_energy_departure_equations(): A.append(coefficients) b.append(value) A, objs = dictionaries2array(A) - values = solve(A, np.array(b).T).T - for obj, value in zip(objs, values): - obj._update_energy_variable(value) + departures = solve(A, np.array(b).T).T + for obj, departure in zip(objs, departures): + if not obj._update_energy_variable: + raise NotImplementedError(f'{obj!r} has no method `_update_energy_variable`') + obj._update_energy_variable(departure) def __enter__(self): units = self.stages @@ -2293,18 +2299,22 @@ def run_phenomena(self): except: for i in path[n+1:]: i.run() last.run() - for i in self.stages: - if (hasattr(i, '_update_equilibrium_variables') - and getattr(i, 'phases', None) == ('g', 'l')): - i._update_equilibrium_variables() - if hasattr(i, '_update_reaction_conversion'): - i._update_reaction_conversion() try: with self.stage_configuration(aggregated=False) as conf: + for i in conf.stages: + if hasattr(i, '_update_equilibrium_variables'): + i._update_equilibrium_variables() + if hasattr(i, '_update_reaction_conversion'): + i._update_reaction_conversion() conf.solve_energy_departures() conf.solve_material_flows() except: with self.stage_configuration(aggregated=True) as conf: + for i in conf.stages: + if hasattr(i, '_update_equilibrium_variables'): + i._update_equilibrium_variables() + if hasattr(i, '_update_reaction_conversion'): + i._update_reaction_conversion() conf.solve_energy_departures() conf.solve_material_flows() diff --git a/biosteam/units/compressor.py b/biosteam/units/compressor.py index d4eeab0c..39e24d53 100644 --- a/biosteam/units/compressor.py +++ b/biosteam/units/compressor.py @@ -600,8 +600,8 @@ def _create_material_balance_equations(self): return equations - def _update_energy_variable(self, variable, value): - self.outs[0].T += value + def _update_energy_variable(self, departure): + self.outs[0].T += departure class PolytropicCompressor(Compressor, new_graphics=False): diff --git a/biosteam/units/distillation.py b/biosteam/units/distillation.py index 2fbe50d8..7b83d524 100644 --- a/biosteam/units/distillation.py +++ b/biosteam/units/distillation.py @@ -1514,6 +1514,12 @@ def _create_material_balance_equations(self): ) return equations + def _get_energy_departure_coefficient(self, stream): + return None + + def _create_energy_departure_equations(self): + return [] + def _update_equilibrium_variables(self): outs = top, bottom = self.outs data = [i.get_data() for i in outs] @@ -1940,6 +1946,8 @@ def _recompute_distillate_recoveries(self, distillate_recoveries): self._distillate_recoveries = distillate_recoveries return distillate_recoveries + _get_energy_departure_coefficient = BinaryDistillation._get_energy_departure_coefficient + _create_energy_departure_equations = BinaryDistillation._create_energy_departure_equations _create_material_balance_equations = BinaryDistillation._create_material_balance_equations _update_equilibrium_variables = BinaryDistillation._update_equilibrium_variables diff --git a/biosteam/units/mixing.py b/biosteam/units/mixing.py index 9104466f..d27bd153 100644 --- a/biosteam/units/mixing.py +++ b/biosteam/units/mixing.py @@ -176,12 +176,12 @@ def _create_material_balance_equations(self): ) return equations - def _update_energy_variable(self, value): + def _update_energy_variable(self, departure): phases = self.phases if phases == ('g', 'l'): - self._B += value + self._B += departure else: - self.outs[0].T += value + self.outs[0].T += departure class SteamMixer(Unit): diff --git a/biosteam/units/splitting.py b/biosteam/units/splitting.py index be89d75e..8b53c00c 100644 --- a/biosteam/units/splitting.py +++ b/biosteam/units/splitting.py @@ -342,7 +342,7 @@ def _create_energy_departure_equations(self, temperature_only=False): self.ins[0]._update_energy_departure_coefficient(coeff) return [(coeff, self.H_in - self.H_out)] - def _update_energy_variable(self, value): - self.T = T = self.T + value + def _update_energy_variable(self, departure): + self.T = T = self.T + departure for i in self.outs: i.T = T diff --git a/biosteam/units/stage.py b/biosteam/units/stage.py index 44803a70..14d468f5 100644 --- a/biosteam/units/stage.py +++ b/biosteam/units/stage.py @@ -102,6 +102,7 @@ def _get_energy_departure_coefficient(self, stream): if self.T is None: return (self, -stream.C) def _create_energy_departure_equations(self): + if self.T is not None: return [] # Ll: C1dT1 - Ce2*dT2 - Cr0*dT0 - hv2*L2*dB2 = Q1 - H_out + H_in # gl: hV1*L1*dB1 - hv2*L2*dB2 - Ce2*dT2 - Cr0*dT0 = Q1 + H_in - H_out outlet = self.outs[0] @@ -129,6 +130,10 @@ def _create_material_balance_equations(self): (eq_overall, sum([i.mol for i in fresh_inlets], zeros)) ] + def _update_energy_variable(self, value): + if self.T is None: return + for i in self.outs: i.T += value + class ReactivePhaseStage(bst.Unit): # Does not include VLE _N_outs = _N_ins = 1 @@ -168,7 +173,7 @@ def _create_material_balance_equations(self): # Overall flows eq_overall = {} predetermined_flow = SparseVector.from_dict(sum_sparse_vectors([i.mol for i in fresh_inlets]), size=n) - rhs = predetermined_flow + reaction._dmol + rhs = predetermined_flow + reaction.dmol eq_overall[product] = ones for i in process_inlets: eq_overall[i] = minus_ones else: @@ -176,7 +181,7 @@ def _create_material_balance_equations(self): # Overall flows eq_overall = {} predetermined_flow = SparseVector.from_dict(sum_sparse_vectors([i.mol for i in fresh_inlets]), size=n) - rhs = predetermined_flow + self._dmol + rhs = predetermined_flow + self.dmol rhs[index] = predetermined_flow[index] if reaction.X == 1: eq_overall[product] = ones @@ -219,7 +224,7 @@ def _create_energy_departure_equations(self): return [(coeff, self.Hnet)] def _update_reaction_conversion(self): - self._dmol = 0.99 * (self._dmol + self.reaction.conversion(self.outs[0])) + self.dmol = 0.99 * (self.dmol + self.reaction.conversion(self.outs[0])) # %% Two phases @@ -408,7 +413,7 @@ def _get_energy_departure_coefficient(self, stream): elif self.T_specification is None: coeff = -stream.C else: - coeff = None + return None return (self, coeff) def _create_energy_departure_equations(self): @@ -475,18 +480,18 @@ def _create_material_balance_equations(self): if reaction: # Reactive liquid if isinstance(reaction, (bst.Rxn, bst.RxnI)): # predetermined_flow = SparseVector.from_dict(sum_sparse_vectors([i.mol for i in fresh_inlets]), size=N) - # rhs = predetermined_flow + self.partition._dmol + # rhs = predetermined_flow + self.partition.dmol # for i in self.outs: eq_overall[i] = ones # for i in process_inlets: eq_overall[i] = minus_ones if reaction._rate: predetermined_flow = SparseVector.from_dict(sum_sparse_vectors([i.mol for i in fresh_inlets]), size=N) - rhs = predetermined_flow + self.partition._dmol + rhs = predetermined_flow + self.partition.dmol for i in self.outs: eq_overall[i] = ones for i in process_inlets: eq_overall[i] = minus_ones else: index = reaction._X_index[1] if reaction.phases else reaction._X_index predetermined_flow = SparseVector.from_dict(sum_sparse_vectors([i.mol for i in fresh_inlets]), size=N) - rhs = predetermined_flow + self.partition._dmol + rhs = predetermined_flow + self.partition.dmol rhs[index] = predetermined_flow[index] if reaction.X == 1: for i in self.outs: eq_overall[i] = ones @@ -553,12 +558,12 @@ def _create_material_balance_equations(self): ) return equations - def _update_energy_variable(self, value): + def _update_energy_variable(self, departure): phases = self.phases if phases == ('g', 'l'): - self.B += value + self.B += departure elif phases == ('L', 'l'): - self.T = T = self.T + value + self.T = T = self.T + departure for i in self.outs: i.T = T else: raise RuntimeError('invalid phases') @@ -571,13 +576,14 @@ def _update_equilibrium_variables(self): T = partition.T for i in (partition.outs, self.outs): i.T = T elif phases == ('L', 'l'): - self.partition._run_lle(update=False) + return + self.partition._run_decoupled_Kgamma() else: raise NotImplementedError(f'K for phases {phases} is not yet implemented') def _update_reaction_conversion(self): if self.reaction and self.phases == ('g', 'l'): - self.partition._run_decoupled_reaction() + self.partition._run_decoupled_reaction() # %% @@ -586,7 +592,9 @@ class PhasePartition(Unit): _N_ins = 1 _N_outs = 2 strict_infeasibility_check = False - reactive_relaxation_factor = 0.5 + dmol_relaxation_factor = 0.5 + T_relaxation_factor = 0.5 + K_relaxation_factor = 0.5 def _init(self, phases, partition_data, top_chemical=None, reaction=None): self.partition_data = partition_data @@ -602,7 +610,7 @@ def _init(self, phases, partition_data, top_chemical=None, reaction=None): self.Q = 0. self.B_specification = self.T_specification = None self.B_fallback = 1 - self._dmol = SparseVector.from_size(self.chemicals.size) + self.dmol = SparseVector.from_size(self.chemicals.size) for i, j in zip(self.outs, self.phases): i.phase = j def _get_mixture(self, linked=True): @@ -650,6 +658,33 @@ def _set_arrays(self, IDs, **kwargs): for i, j in kwargs.items(): setattr(self, i, j) self.IDs = IDs + def _set_arrays_with_relaxation(self, IDs, **kwargs): + IDs_last = self.IDs + IDs = tuple(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(): + last = getattr(self, name) + f = getattr(self, name + '_relaxation_factor') + g = 1 - f + if last.size != size: + last = np.ones(size) + setattr(self, name, last) + for i, j in enumerate(index): + last[j] = f * last[j] + g * array[i] + elif IDs_last and len(IDs_last) < len(IDs): + raise RuntimeError('unknown error') + elif IDs_last == IDs: + for i, j in kwargs.items(): + last = getattr(self, i) + f = getattr(self, i + '_relaxation_factor') + g = 1 - f + last[:] = f * last + g * j + else: + for i, j in kwargs.items(): setattr(self, i, j) + self.IDs = IDs + def _get_activity_model(self): chemicals = self.chemicals index = chemicals.get_lle_indices(sum([i.mol for i in self.ins]).nonzero_keys()) @@ -719,7 +754,9 @@ def _run_decoupled_B(self, stacklevel=1): # Flash Rashford-Rice if phi <= 0 or phi >= 1: return self.B = phi / (1 - phi) - def _run_decoupled_KTvle(self, P=None): # Bubble point + def _run_decoupled_KTvle(self, P=None, + T_relaxation_factor=None, + K_relaxation_factor=None): # Bubble point top, bottom = self.outs if self.T_specification: self._run_vle() @@ -736,13 +773,22 @@ def _run_decoupled_KTvle(self, P=None): # Bubble point x[x == 0] = 1. K_new = p.y / p.x IDs = p.IDs - self.T = p.T - self._set_arrays(IDs, K=K_new) + f = self.T_relaxation_factor if T_relaxation_factor is None else T_relaxation_factor + self.T = f * self.T + (1 - f) * p.T + if K_relaxation_factor is None: + self._set_arrays_with_relaxation(IDs, K=K_new) + else: + krf = self.K_relaxation_factor + self.K_relaxation_factor = K_relaxation_factor + self._set_arrays_with_relaxation(IDs, K=K_new) + self.K_relaxation_factor = krf - def _run_decoupled_reaction(self, P=None): + def _run_decoupled_reaction(self, P=None, relaxation_factor=None): top, bottom = self.outs - f = self.reactive_relaxation_factor - self._dmol = f * self._dmol + (1 - f) * self.reaction.conversion(bottom) + f = self.dmol_relaxation_factor if relaxation_factor is None else relaxation_factor + old = self.dmol + new = self.reaction.conversion(bottom) + self.dmol = f * old + (1 - f) * new def _run_lle(self, P=None, update=True, top_chemical=None): if top_chemical is None: top_chemical = self.top_chemical @@ -799,7 +845,7 @@ def _run_vle(self, P=None, update=True): kwargs['T'] = T ms.vle(**kwargs) index = ms.vle._index - if self.reaction: self._dmol = ms.mol[index] - self.feed.mol[index] + if self.reaction: self.dmol = ms.mol[index] - self.feed.mol[index] IDs = ms.chemicals.IDs IDs = tuple([IDs[i] for i in index]) L_mol = ms.imol['l', IDs] @@ -813,9 +859,10 @@ def _run_vle(self, P=None, update=True): V_total = V_mol.sum() if V_total: y_mol = V_mol / V_total + K_new = y_mol / x_mol else: - y_mol = 0 - K_new = y_mol / x_mol + K_new = np.ones(len(index)) * 1e-16 + if B is None: if not L_total: self.B = inf @@ -1349,7 +1396,7 @@ def _get_energy_departure_coefficient(self, stream): assert self.aggregated if stream.phase != 'g': return None vapor, liquid = self.outs - assert stream is vapor + assert stream.imol is vapor.imol if vapor.isempty(): liquid.phase = 'g' coeff = liquid.H @@ -1383,22 +1430,26 @@ def _create_energy_departure_equations(self): def _create_material_balance_equations(self): top, bottom = self.outs - if bottom.isempty(): - B = np.inf - K = 1e16 * np.ones(self.chemicals.size) - elif top.isempty(): - K = np.zeros(self.chemicals.size) - B = 0 - else: - top_mol = top.mol.to_array() - bottom_mol = bottom.mol.to_array() - F_top = top_mol.sum() - F_bottom = bottom_mol.sum() - y = top_mol / F_top - x = bottom_mol / F_bottom - x[x <= 0] = 1e-16 - K = y / x - B = F_top / F_bottom + try: + B = self.B + K = self.K + except: + if bottom.isempty(): + self.B = B = np.inf + self.K = K = 1e16 * np.ones(self.chemicals.size) + elif top.isempty(): + self.K = K = np.zeros(self.chemicals.size) + self.B = B = 0 + else: + top_mol = top.mol.to_array() + bottom_mol = bottom.mol.to_array() + F_top = top_mol.sum() + F_bottom = bottom_mol.sum() + y = top_mol / F_top + x = bottom_mol / F_bottom + x[x <= 0] = 1e-16 + self.K = K = y / x + self.B = B = F_top / F_bottom inlets = self.ins fresh_inlets = [i for i in inlets if i.isfeed() and not i.material_equations] @@ -1417,9 +1468,16 @@ def _create_material_balance_equations(self): del eq_overall[i] else: eq_overall[i] = minus_ones - equations.append( - (eq_overall, sum([i.mol for i in fresh_inlets], zeros)) - ) + if self.stage_reactions: + partitions = self.partitions + flows = [i.mol for i in fresh_inlets] + [partitions[i].dmol for i in self.stage_reactions] + equations.append( + (eq_overall, sum(flows, zeros)) + ) + else: + equations.append( + (eq_overall, sum([i.mol for i in fresh_inlets], zeros)) + ) # Top to bottom flows eq_outs = {} @@ -1436,15 +1494,46 @@ def _create_material_balance_equations(self): return equations def _update_equilibrium_variables(self): - outs = top, bottom = self.outs self.run() + top, bottom = self.outs + if bottom.isempty(): + self.K = 1e16 * np.ones(self.chemicals.size) + elif top.isempty(): + self.K = np.zeros(self.chemicals.size) + else: + top_mol = top.mol.to_array() + bottom_mol = bottom.mol.to_array() + F_top = top_mol.sum() + F_bottom = bottom_mol.sum() + y = top_mol / F_top + x = bottom_mol / F_bottom + x[x <= 0] = 1e-16 + self.K = y / x - def _update_energy_variable(self, value): + def _update_energy_variable(self, departure): phases = self.phases if phases == ('g', 'l'): - self.B += value + if not hasattr(self, 'B'): + top, bottom = self.outs + if bottom.isempty(): + self.B = np.inf + self.K = 1e16 * np.ones(self.chemicals.size) + elif top.isempty(): + self.K = np.zeros(self.chemicals.size) + self.B = 0 + else: + top_mol = top.mol.to_array() + bottom_mol = bottom.mol.to_array() + F_top = top_mol.sum() + F_bottom = bottom_mol.sum() + y = top_mol / F_top + x = bottom_mol / F_bottom + x[x <= 0] = 1e-16 + self.K = y / x + self.B = F_top / F_bottom + self.B += departure elif phases == ('L', 'l'): - for i in self.outs: i.T += value + for i in self.outs: i.T += departure else: raise RuntimeError('invalid phases') @@ -1471,7 +1560,7 @@ def correct_overall_mass_balance(self): stage_reactions = self.stage_reactions if stage_reactions: partitions = self.partitions - inmol += sum([partitions[i]._dmol for i in stage_reactions]) + inmol += sum([partitions[i].dmol for i in stage_reactions]) try: factor = inmol / outmol except: @@ -1486,7 +1575,7 @@ def material_errors(self): for stage in stages: errors.append( sum([i.mol for i in stage.ins], - -sum([i.mol for i in stage.outs], -stage.partition._dmol)) + -sum([i.mol for i in stage.outs], -stage.partition.dmol)) ) return pd.DataFrame(errors, columns=IDs) @@ -1496,7 +1585,7 @@ def _feed_flows_and_conversion(self): index = self._update_index for i in self.stage_reactions: p = partition[i] - dmol = p._dmol + dmol = p.dmol feed_mol = p.feed.mol pseudo_feed = feed_mol + dmol mask = pseudo_feed < 0 @@ -2111,7 +2200,10 @@ def interpolate_missing_variables(self): B = partition.B T = partition.T K = partition.K - if B is None or K is None or K.size != N_chemicals: continue + try: + if B is None or K is None or K.size != N_chemicals: continue + except: + breakpoint() index.append(i) Bs.append(B) Ks.append(K) diff --git a/tests/test_equilibrium_stages.py b/tests/test_equilibrium_stages.py index 5a67e11f..00b0002d 100644 --- a/tests/test_equilibrium_stages.py +++ b/tests/test_equilibrium_stages.py @@ -194,11 +194,11 @@ def test_distillation(): distillation.simulate() flows = [ [0.0, 0.0, 0.0], - [0.11207831812626445, 4.350506202064179, 22.358675959386957], - [22.055921681873734, 0.14389379793582036, 87.15232404061304], + [0.1120748575014966, 4.350505208368632, 22.358805264239038], + [22.055925142498495, 0.14389479163136612, 87.1521947357609], ] for i, j in zip(distillation.outs, flows): - assert_allclose(i.mol, j, rtol=1e-6, atol=1e-3) + assert_allclose(i.mol, j, rtol=1e-5, atol=1e-3) if __name__ == '__main__': diff --git a/tests/test_phenomena_oriented_simulation.py b/tests/test_phenomena_oriented_simulation.py index d7b97d6f..0e19faa3 100644 --- a/tests/test_phenomena_oriented_simulation.py +++ b/tests/test_phenomena_oriented_simulation.py @@ -207,7 +207,7 @@ def fresh_solvent_flow_rate(): for s_sm, s_dp in zip(sm.streams, po.streams): actual = s_sm.mol value = s_dp.mol - assert_allclose(actual, value, rtol=1e-1, atol=1e-3) + assert_allclose(actual, value, rtol=1e-5, atol=1e-3) if __name__ == '__main__': test_trivial_lle_case()