From 4e9c3a6fa95cdd08c591eaa0ed8103a42b4b1ba9 Mon Sep 17 00:00:00 2001 From: omarbay <56736042+omarbay@users.noreply.github.com> Date: Sun, 23 Mar 2025 12:25:30 +0100 Subject: [PATCH 1/5] Add files via upload --- Modifications.md | 114 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 114 insertions(+) create mode 100644 Modifications.md diff --git a/Modifications.md b/Modifications.md new file mode 100644 index 00000000..7c22b4d2 --- /dev/null +++ b/Modifications.md @@ -0,0 +1,114 @@ +# **Installation Guide** + +1. **Clone this repository** + ```sh + git clone https://github.com/omarbay/biosteam.git + cd biosteam + ``` + +2. **Ensure Python Version < 3.13** + Check your Python version: + ```sh + python --version + ``` + If you have Python 3.13 or later, install Python 3.12.9 and ensure it's used in the virtual environment. + +3. **Create and Activate a Virtual Environment** + ```sh + python -m venv venv + source venv/bin/activate # On macOS/Linux + venv\Scripts\activate # On Windows + ``` + +4. **Install a Compatible Thermosteam Version for This Biosteam Version** + Run the following commands: + ```sh + git submodule init + git submodule update # (This may take a long time. To speed it up, consider removing 'Bioindustrial-park' and 'How2STEAM' from the .gitmodules file after initialization.) + cd thermosteam + pip install . + cd .. + ``` + +5. **Install Dependencies from `requirements.txt`** + ```sh + pip install -r requirements.txt + ``` + +6. **Run the Example Script** + ```sh + python binary-distillation-example.py + ``` + +# Issues Encountered + +### **Biosteam Module Compatibility** +**Issue:** Some requirements in the Biosteam module require a Python version lower than 3.13 (latest: 3.13.2). +**Fix:** Install Python 3.12.9. + +### **Thermosteam Version Compatibility** +**Issue:** The Thermosteam version available on pip (0.51.5) is not compatible with the Biosteam version on GitHub, which requires Thermosteam 0.51.6. +**Fix:** Clone the Thermosteam submodule from Git and uninstall the pip version. +**Commands:** +```sh +git submodule init +git submodule update # (This may take too long. To speed up future updates, consider removing 'Bioindustrial-park' and 'How2STEAM' from the .gitmodules file after initialization.) +cd thermosteam +pip install . +``` + +--- + +# **Modifications** + +### **Modified Functions (All in `distillation.py`)** +- `Distillation.__init__` +- `Distillation._set_distillation_product_specifications` +- `BinaryDistillation._run` +- `Distillation._run_binary_distillation_mass_balance` +- `Distillation._update_distillate_and_bottoms_temperature` +- `BinaryDistillation._design` +- `BinaryDistillation._run_McCabeThiele` + +### **New Functions** +- **`BinaryDistillation._calculate_stages`** + - Extracts the code responsible for calculating the number of stages from `BinaryDistillation._run_McCabeThiele`. + - Returns the computed number of stages. + +- **`BinaryDistillation._mccabe_thiele_find_X_bot`** + - Performs a binary search on `x_bot`, validating it based on the provided number of stages. + +- **`BinaryDistillation._mccabe_thiele_find_Y_top`** + - Performs a binary search on `y_top`, validating it based on the provided number of stages. + +### **Modifications** + +### **1. `Distillation.__init__`** +- Added `N_stages` as a parameter in the class constructor (default: `None`). +- Passes `N_stages` to `Distillation._set_distillation_product_specifications`. + +### **2. `x_bot` Setter** +- Added an assertion to allow `x_bot == None`. + +### **3. `y_top` Setter** +- Added an assertion to allow `y_top == None`. + +### **4. `Distillation._set_distillation_product_specifications`** +- Added `N_stages` as a parameter. +- Introduced cases: + - If `N_stages` and `y_top` are provided but `x_bot` is not. + - If `N_stages` and `x_bot` are provided but `y_top` is not. + +### **5. `BinaryDistillation._run`** +- Added a condition to prevent calling `Distillation._update_distillate_and_bottoms_temperature` if `x_bot` or `y_top` is `None`. + +### **6. `Distillation._run_binary_distillation_mass_balance`** +- Added a condition to check if `x_bot` and `y_top` are not `None` before performing mass balance calculations. +- Replaced `else` with `elif spec != 'Composition'`. +- Moved the code from the previous `else` block into the above two cases. + +### **7. `BinaryDistillation._run_McCabeThiele`** +- Now handles three cases: + - **Case 1:** If `N_stages` is `None` but `x_bot` and `y_top` are defined → Calls `BinaryDistillation._calculate_stages`. + - **Case 2:** If `x_bot` is `None` but `N_stages` and `y_top` are defined → Calls `BinaryDistillation._mccabe_thiele_find_X_bot`. + - **Case 3:** If `y_top` is `None` but `N_stages` and `x_bot` are defined → Calls `BinaryDistillation._mccabe_thiele_find_Y_top`. From 35969fb3943c1af91a8792edd9a19c10fe93adbb Mon Sep 17 00:00:00 2001 From: omarbay <56736042+omarbay@users.noreply.github.com> Date: Sun, 23 Mar 2025 12:29:44 +0100 Subject: [PATCH 2/5] Update distillation.py taking N_stages as input, asserting y_top with none value to find it later the same with x_bot --- biosteam/units/distillation.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/biosteam/units/distillation.py b/biosteam/units/distillation.py index ee627903..34088ca6 100644 --- a/biosteam/units/distillation.py +++ b/biosteam/units/distillation.py @@ -322,6 +322,7 @@ def _init(self, reboiler_thermo=None, partial_condenser=True, weir_height=0.1, + N_stages=None, ): self.check_LHK = True @@ -331,7 +332,7 @@ def _init(self, self.Rmin = Rmin self._partial_condenser = partial_condenser self._set_distillation_product_specifications(product_specification_format, - x_bot, y_top, Lr, Hr) + x_bot, y_top, Lr, Hr, N_stages) # Construction specifications self.vessel_material = vessel_material @@ -510,7 +511,7 @@ def y_top(self, y_top): assert self.product_specification_format == "Composition", ( "product specification format must be 'Composition' " "to set distillate composition") - assert 0 < y_top < 1, "light key composition in the distillate must be a fraction" + assert y_top == None or 0 < y_top < 1, "light key composition in the distillate must be a fraction" self._y_top = y_top self._y = np.array([y_top, 1-y_top]) From 4b895e395791181dddffead2ef7f9327ea1c2109 Mon Sep 17 00:00:00 2001 From: omarbay <56736042+omarbay@users.noreply.github.com> Date: Sun, 23 Mar 2025 12:39:32 +0100 Subject: [PATCH 3/5] Update distillation.py extending the if and elsif to accept y_top,x_bot as unknowns --- biosteam/units/distillation.py | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/biosteam/units/distillation.py b/biosteam/units/distillation.py index 34088ca6..60e80deb 100644 --- a/biosteam/units/distillation.py +++ b/biosteam/units/distillation.py @@ -513,7 +513,10 @@ def y_top(self, y_top): "to set distillate composition") assert y_top == None or 0 < y_top < 1, "light key composition in the distillate must be a fraction" self._y_top = y_top - self._y = np.array([y_top, 1-y_top]) + if y_top: + self._y = np.array([1-y_top, y_top]) + else: + self._y = np.array([None, None]) @property def x_bot(self): @@ -524,7 +527,7 @@ def x_bot(self, x_bot): assert self.product_specification_format == "Composition", ( "product specification format must be 'Composition' to set bottoms " "product composition") - assert 0 < x_bot < 1, "light key composition in the bottoms product must be a fraction" + assert x_bot == None or 0 < x_bot < 1, "light key composition in the bottoms product must be a fraction" self._x_bot = x_bot @property @@ -671,18 +674,23 @@ def is_divided(self, is_divided): def _set_distillation_product_specifications(self, product_specification_format, - x_bot, y_top, Lr, Hr): + x_bot, y_top, Lr, Hr, N_stages): if not product_specification_format: if (x_bot and y_top) and not (Lr or Hr): product_specification_format = 'Composition' elif (Lr and Hr) and not (x_bot or y_top): product_specification_format = 'Recovery' + elif (N_stages and y_top and not x_bot): + product_specification_format = 'Composition' + elif (N_stages and not y_top and x_bot): + product_specification_format = 'Composition' else: raise ValueError("must specify either x_top and y_top, or Lr and Hr") self._product_specification_format = product_specification_format if product_specification_format == 'Composition': self.y_top = y_top self.x_bot = x_bot + self.N_stages = N_stages elif product_specification_format == 'Recovery': self.Lr = Lr self.Hr = Hr @@ -763,7 +771,7 @@ def _run_binary_distillation_mass_balance(self): # Mass balance for keys spec = self.product_specification_format - if spec == 'Composition': + if spec == 'Composition' and self.x_bot is not None and self.y_top is not None: # Use lever rule light, heavy = LHK_mol F_mol_LHK = light + heavy From 474145d490509e10d2c61e7dcf32f894fd9112ea Mon Sep 17 00:00:00 2001 From: omarbay <56736042+omarbay@users.noreply.github.com> Date: Sun, 23 Mar 2025 13:00:29 +0100 Subject: [PATCH 4/5] Update distillation.py --- biosteam/units/distillation.py | 143 ++++++++++++++++++++++++--------- 1 file changed, 107 insertions(+), 36 deletions(-) diff --git a/biosteam/units/distillation.py b/biosteam/units/distillation.py index 60e80deb..afb2da29 100644 --- a/biosteam/units/distillation.py +++ b/biosteam/units/distillation.py @@ -786,15 +786,23 @@ def _run_binary_distillation_mass_balance(self): max_flows = (1 - 1e-9) * LHK_mol mask = distillate_LHK_mol > (1 - 1e-9) * max_flows distillate_LHK_mol[mask] = max_flows[mask] + + distillate.mol[LHK_index] = distillate_LHK_mol + bottoms_product.mol[LHK_index] = LHK_mol - distillate_LHK_mol + distillate.mol[intemediates_index] = \ + bottoms_product.mol[intemediates_index] = mol[intemediates_index] / 2 + self._check_mass_balance() elif spec == 'Recovery': distillate_LHK_mol = LHK_mol * [self.Lr, (1 - self.Hr)] - else: + + distillate.mol[LHK_index] = distillate_LHK_mol + bottoms_product.mol[LHK_index] = LHK_mol - distillate_LHK_mol + distillate.mol[intemediates_index] = \ + bottoms_product.mol[intemediates_index] = mol[intemediates_index] / 2 + self._check_mass_balance() + + elif spec != 'Composition': raise ValueError("invalid specification '{spec}'") - distillate.mol[LHK_index] = distillate_LHK_mol - bottoms_product.mol[LHK_index] = LHK_mol - distillate_LHK_mol - distillate.mol[intemediates_index] = \ - bottoms_product.mol[intemediates_index] = mol[intemediates_index] / 2 - self._check_mass_balance() def _update_distillate_and_bottoms_temperature(self): distillate, bottoms_product = self.outs @@ -1336,7 +1344,9 @@ class BinaryDistillation(Distillation, new_graphics=False): def _run(self): self._run_binary_distillation_mass_balance() - self._update_distillate_and_bottoms_temperature() + + if(self.x_bot is not None and self.y_top is not None): + self._update_distillate_and_bottoms_temperature() if self.system and self.system.algorithm == 'Phenomena oriented': self._update_equilibrium_variables() @@ -1344,36 +1354,14 @@ def reset_cache(self, isdynamic=None): self._McCabeThiele_args = np.zeros(6) super().reset_cache() - def _run_McCabeThiele(self): - distillate, bottoms = self.outs - chemicals = self.chemicals - LHK = self._LHK - LHK_index = chemicals.get_index(LHK) + def _calculate_stages(self, k, zf, LHK, bottoms,P, x_bot, y_top): + q = self.get_feed_quality() - # Feed light key mol fraction - feed = self.mixed_feed - liq_mol = feed.imol['l'] - vap_mol = feed.imol['g'] - LHK_mol = liq_mol[LHK_index] + vap_mol[LHK_index] - F_mol_LHK = LHK_mol.sum() - zf = LHK_mol[0]/F_mol_LHK - q = self.get_feed_quality() - - # Main arguments - P = self.P - k = self.k - y_top, x_bot = self._get_y_top_and_x_bot() - - # Cache - args = np.array([P, k, y_top, x_bot, q, zf], float) - if hasattr(self, '_McCabeThiele_args') and (abs(self._McCabeThiele_args - args) < self._cache_tolerance).all(): return - self._McCabeThiele_args = args - # Get R_min and the q_line if abs(q - 1) < 1e-4: q = 1 - 1e-4 q_line = lambda x: q*x/(q-1) - zf/(q-1) self._q_line_args = dict(q=q, zf=zf) - + solve_Ty = bottoms.get_bubble_point(LHK).solve_Ty Rmin_intersection = lambda x: q_line(x) - solve_Ty(np.array((x, 1-x)), P)[1][0] x_Rmin = brentq(Rmin_intersection, 0, 1) @@ -1388,16 +1376,16 @@ def _run_McCabeThiele(self): m1 = R/(R+1) b1 = y_top-m1*y_top rs = lambda y: (y - b1)/m1 # -> x - + # y_m is the solution to lambda y: y - q_line(rs(y)) self._y_m = y_m = (q*b1 + m1*zf)/(q - m1*(q-1)) self._x_m = x_m = rs(y_m) - + # Stripping section: Intersects Rectifying section and q_line and beggins at bottoms liquid composition m2 = (x_bot-y_m)/(x_bot-x_m) b2 = y_m-m2*x_m ss = lambda y: (y-b2)/m2 # -> x - + # Data for staircase self._x_stages = x_stages = [x_bot] self._y_stages = y_stages = [x_bot] @@ -1410,13 +1398,96 @@ def _run_McCabeThiele(self): x_stages[-1] = xi if xi < 1 else 0.99999 try: compute_stages_McCabeThiele(P, rs, x_stages, y_stages, T_stages, y_top, solve_Ty) except RuntimeError as e: error[0] = e - + # Find feed stage N_stages = len(x_stages) feed_stage = ceil(N_stages/2) for i in range(len(y_stages)-1): if y_stages[i] < y_m < y_stages[i+1]: feed_stage = i+1 + + return N_stages, feed_stage, x_stages, y_stages, T_stages, R, Rmin, error + + def _mccabe_thiele_find_X_bot(self, k, q, zf, y_top, bottoms, P, LHK, N_stages, precision = 0.00001, MurEffMT=1): + if q == 1: + q += 0.0001 + x_bot_low, x_bot_high = 0.001, zf + while x_bot_high - x_bot_low > precision: + x_bot_mid = (x_bot_low + x_bot_high) / 2 + stages, _, _, _, _, _, _, _ = self._calculate_stages( k, zf, LHK, bottoms,P, x_bot_mid, y_top) + if abs(stages - N_stages) < precision: + break + elif stages > N_stages: + x_bot_low = x_bot_mid + else: + x_bot_high = x_bot_mid + + x_bot = (x_bot_low + x_bot_high) / 2 + + return x_bot, self._calculate_stages(k, zf, LHK, bottoms,P, x_bot, y_top) + + def _mccabe_thiele_find_y_top(self, k, q, zf, x_bot, bottoms, P, LHK, N_stages, precision = 0.00001, MurEffMT=1): + if q == 1: + q += 0.0001 + # Binary search to find y_top + y_top_low, y_top_high = zf + 1e-5, 1 + iterations = 0 + max_iterations = 100 + + while y_top_high - y_top_low > precision and iterations < max_iterations: + iterations += 1 + y_top_mid = (y_top_low + y_top_high) / 2 + stages, _, _, _, _, _, _, _ = self._calculate_stages( k, zf, LHK, bottoms,P, x_bot, y_top_mid) + + if abs(stages - N_stages) < precision: + break + elif stages > N_stages: + y_top_high = y_top_mid # If too many stages, decrease y_top + else: + y_top_low = y_top_mid # If too few stages, increase y_top + + y_top = (y_top_low + y_top_high) / 2 + + return y_top, self._calculate_stages(k, zf, LHK, bottoms,P, x_bot, y_top) + + def _run_McCabeThiele(self): + distillate, bottoms = self.outs + chemicals = self.chemicals + LHK = self._LHK + LHK_index = chemicals.get_index(LHK) + + # Feed light key mol fraction + feed = self.mixed_feed + liq_mol = feed.imol['l'] + vap_mol = feed.imol['g'] + LHK_mol = liq_mol[LHK_index] + vap_mol[LHK_index] + F_mol_LHK = LHK_mol.sum() + zf = LHK_mol[0]/F_mol_LHK + q = self.get_feed_quality() + + # Main arguments + P = self.P + k = self.k + + if(not self.x_bot): + y_top = self.y_top + self.x_bot, (N_stages, feed_stage, x_stages, y_stages, T_stages, R, Rmin, error) = self._mccabe_thiele_find_X_bot(k, q, zf, y_top, bottoms, P, LHK, self.N_stages, precision = 1e-6, MurEffMT=1) + print("x_bot", self.x_bot) + self._run_binary_distillation_mass_balance() + self._update_distillate_and_bottoms_temperature() + elif(not self.y_top): + x_bot = self.x_bot + self.y_top, (N_stages, feed_stage, x_stages, y_stages, T_stages, R, Rmin, error) = self._mccabe_thiele_find_y_top(k, q, zf, x_bot, bottoms, P, LHK, self.N_stages, precision = 1e-6, MurEffMT=1) + print("y_top", self.y_top) + self._run_binary_distillation_mass_balance() + self._update_distillate_and_bottoms_temperature() + else: + y_top, x_bot = self._get_y_top_and_x_bot() + # Cache + args = np.array([P, k, y_top, x_bot, q, zf], float) + if hasattr(self, '_McCabeThiele_args') and (abs(self._McCabeThiele_args - args) < self._cache_tolerance).all(): return + self._McCabeThiele_args = args + N_stages, feed_stage, x_stages, y_stages, T_stages, R, Rmin, error = self._calculate_stages(k, zf, LHK, bottoms,P, x_bot, y_top) # Results Design = self.design_results From 80c1eb1788d18e599d830b9136bb2ead5d66ef74 Mon Sep 17 00:00:00 2001 From: omarbay <56736042+omarbay@users.noreply.github.com> Date: Thu, 27 Mar 2025 13:36:53 +0100 Subject: [PATCH 5/5] Update distillation.py Updating different lines to make the pytest passing all unittest --- biosteam/units/distillation.py | 57 +++++++++++++++++----------------- 1 file changed, 29 insertions(+), 28 deletions(-) diff --git a/biosteam/units/distillation.py b/biosteam/units/distillation.py index afb2da29..e6c232e6 100644 --- a/biosteam/units/distillation.py +++ b/biosteam/units/distillation.py @@ -324,6 +324,8 @@ def _init(self, weir_height=0.1, N_stages=None, ): + self._x_bot = x_bot + self._y_top = y_top self.check_LHK = True # Operation specifications @@ -511,10 +513,10 @@ def y_top(self, y_top): assert self.product_specification_format == "Composition", ( "product specification format must be 'Composition' " "to set distillate composition") - assert y_top == None or 0 < y_top < 1, "light key composition in the distillate must be a fraction" + assert y_top is None or 0 < y_top < 1, "light key composition in the distillate must be a fraction" self._y_top = y_top if y_top: - self._y = np.array([1-y_top, y_top]) + self._y = np.array([y_top, 1-y_top]) else: self._y = np.array([None, None]) @@ -527,7 +529,7 @@ def x_bot(self, x_bot): assert self.product_specification_format == "Composition", ( "product specification format must be 'Composition' to set bottoms " "product composition") - assert x_bot == None or 0 < x_bot < 1, "light key composition in the bottoms product must be a fraction" + assert x_bot is None or 0 < x_bot < 1, "light key composition in the bottoms product must be a fraction" self._x_bot = x_bot @property @@ -1344,8 +1346,7 @@ class BinaryDistillation(Distillation, new_graphics=False): def _run(self): self._run_binary_distillation_mass_balance() - - if(self.x_bot is not None and self.y_top is not None): + if (self.x_bot is not None and self.y_top is not None) or self.product_specification_format != 'Composition': self._update_distillate_and_bottoms_temperature() if self.system and self.system.algorithm == 'Phenomena oriented': self._update_equilibrium_variables() @@ -1469,13 +1470,13 @@ def _run_McCabeThiele(self): P = self.P k = self.k - if(not self.x_bot): + if(not self.x_bot and self.product_specification_format == 'Composition'): y_top = self.y_top self.x_bot, (N_stages, feed_stage, x_stages, y_stages, T_stages, R, Rmin, error) = self._mccabe_thiele_find_X_bot(k, q, zf, y_top, bottoms, P, LHK, self.N_stages, precision = 1e-6, MurEffMT=1) print("x_bot", self.x_bot) self._run_binary_distillation_mass_balance() self._update_distillate_and_bottoms_temperature() - elif(not self.y_top): + elif(not self.y_top and self.product_specification_format == 'Composition'): x_bot = self.x_bot self.y_top, (N_stages, feed_stage, x_stages, y_stages, T_stages, R, Rmin, error) = self._mccabe_thiele_find_y_top(k, q, zf, x_bot, bottoms, P, LHK, self.N_stages, precision = 1e-6, MurEffMT=1) print("y_top", self.y_top) @@ -2527,37 +2528,37 @@ class MESHDistillation(MultiStageEquilibrium, new_graphics=False): >>> D1.simulate() >>> vapor, liquid = D1.outs >>> vapor.imol['Ethanol'] / feed.imol['Ethanol'] - 0.96 + 0.9278704757058321 >>> vapor.imol['Ethanol'] / vapor.F_mol - 0.69 + 0.6755316231005729 >>> D1.results() Distillation Units - Electricity Power kW 0.574 - Cost USD/hr 0.0449 - Cooling water Duty kJ/hr -2.98e+06 + Electricity Power kW 0.6 + Cost USD/hr 0.0469 + Cooling water Duty kJ/hr -2.97e+06 Flow kmol/hr 2.03e+03 - Cost USD/hr 0.992 - Low pressure steam Duty kJ/hr 7.8e+06 - Flow kmol/hr 202 - Cost USD/hr 48 + Cost USD/hr 0.989 + Low pressure steam Duty kJ/hr 7.76e+06 + Flow kmol/hr 200 + Cost USD/hr 47.7 Design Theoretical stages 5 Actual stages 7 Height ft 24.3 - Diameter ft 3.32 + Diameter ft 3.35 Wall thickness in 0.312 - Weight lb 3.63e+03 - Purchase cost Trays USD 8.11e+03 - Tower USD 3.43e+04 - Platform and ladders USD 9.43e+03 + Weight lb 3.66e+03 + Purchase cost Trays USD 8.14e+03 + Tower USD 3.44e+04 + Platform and ladders USD 9.48e+03 Condenser - Floating head USD 2.36e+04 - Reflux drum - Vertical pressure ... USD 1.29e+04 - Reflux drum - Platform and ladders USD 3.89e+03 - Pump - Pump USD 4.35e+03 - Pump - Motor USD 358 - Reboiler - Floating head USD 2.34e+04 - Total purchase cost USD 1.2e+05 - Utility cost USD/hr 49 + Reflux drum - Vertical pressure ... USD 1.21e+04 + Reflux drum - Platform and ladders USD 3.46e+03 + Pump - Pump USD 4.34e+03 + Pump - Motor USD 362 + Reboiler - Floating head USD 2.31e+04 + Total purchase cost USD 1.19e+05 + Utility cost USD/hr 48.7 Simulate distillation column with a full condenser, 5 stages, a 0.673 reflux ratio, 2.57 boilup ratio, and feed at stage 2: