From 321bba771b8d85db9e3af6d372c8ed6f0e56b45c Mon Sep 17 00:00:00 2001 From: daklauss Date: Wed, 31 Jul 2024 09:32:12 +0200 Subject: [PATCH] implemented componentsystem that inherits from cadet process and fixes the small bug with the molecular weight, added molecular volume first implementation of the def unit_operation and residual equations still both still need testing and unit_operation def doesn't have rejection model yet --- CADETPythonSimulator/componentsystem.py | 265 ++++++++++++++++++++++++ CADETPythonSimulator/residual.py | 33 +-- CADETPythonSimulator/system_solver.py | 1 + CADETPythonSimulator/unit_operation.py | 100 ++++++--- pyproject.toml | 6 +- tests/test_residual.py | 39 +++- tests/test_unit_operation.py | 22 +- 7 files changed, 394 insertions(+), 72 deletions(-) create mode 100644 CADETPythonSimulator/componentsystem.py diff --git a/CADETPythonSimulator/componentsystem.py b/CADETPythonSimulator/componentsystem.py new file mode 100644 index 0000000..d790180 --- /dev/null +++ b/CADETPythonSimulator/componentsystem.py @@ -0,0 +1,265 @@ +from CADETProcess.processModel import ComponentSystem, Component, Species +from CADETProcess.dataStructure import UnsignedFloat, String, Integer +from CADETProcess.dataStructure import Structure + +from CADETPythonSimulator.exception import CADETPythonSimError +from functools import wraps + +class CPSSpecies(Structure): + """Species class. + + Represent a species in a chemical system. + Same as in cadet Process but with added + density and Volume + + Attributes + ---------- + name : str + The name of the species. + charge : int, optional + The charge of the species. Default is 0. + molecular_weight : float + The molecular weight of the species. + density : float + Density of the species. + molecular_volume : float + The molecular volume of the species + + """ + name = String() + charge = Integer(default=0) + molecular_weight = UnsignedFloat() + density = UnsignedFloat() + molecular_volume = UnsignedFloat() + +class CPSComponent(Component): + """Information about single component. + Inherits from CadetProcess Component + Same function but with fixed molecular weight and added densities and volume + + A component can contain subspecies (e.g. differently charged variants). + + Attributes + ---------- + name : String + Name of the component. + species : list + List of Subspecies. + n_species : int + Number of Subspecies. + label : list + Name of component (including species). + charge : list + Charge of component (including species). + molecular_weight : list + Molecular weight of component (including species). + density : list + density of component (including species). + molecular_volume : list + Molecular volume of component (including species). + + See Also + -------- + Species + ComponentSystem + + """ + def __init__(self, + name=None, + species=None, + charge=None, + molecular_weight=None, + density=None, + molecular_volume=None): + + self.name = name + self._species = [] + + if species is None: + self.add_species(name, charge, molecular_weight, density, molecular_volume) + elif isinstance(species, str): + self.add_species(species, + charge, molecular_weight, density, molecular_volume) + elif isinstance(species, list): + if charge is None: + charge = len(species) * [None] + if molecular_weight is None: + molecular_weight = len(species) * [None] + if density is None: + density = len(species) * [None] + if molecular_volume is None: + molecular_volume = len(species) * [None] + for i, spec in enumerate(species): + self.add_species(spec, + charge[i], molecular_weight[i], density[i], molecular_volume[i]) + else: + raise CADETPythonSimError("Could not determine number of species") + + def add_species(self, species, *args, **kwargs): + if not isinstance(species, CPSSpecies): + species = CPSSpecies(species, *args, **kwargs) + self._species.append(species) + + @property + def molecular_volume(self): + """list of float or None: The molecular volume of the subspecies.""" + return [spec.molecular_volume for spec in self.species] + + @property + def density(self): + """list of float or None: The density of the subspecies.""" + return [spec.density for spec in self.species] + + @property + def molecular_weight(self): + """list of float or None: The molecular weights of the subspecies.""" + return [spec.molecular_weight for spec in self.species] + +class CPSComponentSystem(ComponentSystem): + """Information about components in system. Inherits from Component System. Adds + molecular Volume to the Component System. + + A component can contain subspecies (e.g. differently charged variants). + + Attributes + ---------- + name : String + Name of the component system. + components : list + List of individual components. + n_species : int + Number of Subspecies. + n_comp : int + Number of all component species. + n_components : int + Number of components. + indices : dict + Component indices. + names : list + Names of all components. + species : list + Names of all component species. + charge : list + Charges of all components species. + molecular_weight : list + Molecular weights of all component species. + molecular_volume : list + Molecular volume of all component species. + + See Also + -------- + Species + Component + + """ + + def __init__( + self, + components=None, + name=None, + charges=None, + molecular_weights=None, + densities=None, + molecular_volume=None + ): + """Initialize the ComponentSystem object. + + Parameters + ---------- + components : int, list, None + The number of components or the list of components to be added. + If None, no components are added. + name : str, None + The name of the ComponentSystem. + charges : list, None + The charges of each component. + molecular_weights : list, None + The molecular weights of each component. + densities : list, None + Densities of each component + molecular_volume : list, None + The molecular volume of each component. + + Raises + ------ + CADETProcessError + If the `components` argument is neither an int nor a list. + + """ + + self.name = name + + self._components = [] + + if components is None: + return + + if isinstance(components, int): + n_comp = components + components = [str(i) for i in range(n_comp)] + elif isinstance(components, list): + n_comp = len(components) + else: + raise CADETPythonSimError("Could not determine number of components") + + if charges is None: + charges = n_comp * [None] + if molecular_weights is None: + molecular_weights = n_comp * [None] + if densities is None: + densities = n_comp * [None] + if molecular_volume is None: + molecular_volume = n_comp * [None] + + for i, comp in enumerate(components): + self.add_component( + comp, + charge=charges[i], + molecular_weight=molecular_weights[i], + density=densities[i], + molecular_volume=molecular_volume[i] + ) + + @wraps(CPSComponent.__init__) + def add_component(self, component, *args, **kwargs): + """ + Add a component to the system. + + Parameters + ---------- + component : {str, Component} + The class of the component to be added. + *args : list + The positional arguments to be passed to the component class's constructor. + **kwargs : dict + The keyword arguments to be passed to the component class's constructor. + + """ + if not isinstance(component, CPSComponent): + component = CPSComponent(component, *args, **kwargs) + + if component.name in self.names: + raise CADETPythonSimError( + f"Component '{component.name}' " + "already exists in ComponentSystem." + ) + + self._components.append(component) + + @property + def molecular_volumes(self): + """list: List of species molecular volumes.""" + molecular_volumes = [] + for comp in self.components: + molecular_volumes += comp.molecular_volume + + return molecular_volumes + + @property + def densities(self): + """list: List of species densities.""" + densities = [] + for comp in self.components: + densities += comp.density + + return densities diff --git a/CADETPythonSimulator/residual.py b/CADETPythonSimulator/residual.py index da87b03..8fc2318 100644 --- a/CADETPythonSimulator/residual.py +++ b/CADETPythonSimulator/residual.py @@ -23,10 +23,10 @@ def calculate_residual_volume_cstr( Volume leaving the Unit Returns ------- - float + float Residual of the Flow equation of the CSTR with dimensions like the inpu """ - + if V < 0: raise CADETPythonSimError("V can't be less then zero") @@ -43,7 +43,7 @@ def calculate_residual_concentration_cstr( ) -> np.ndarray : """ Calculates the residual equations of the concentration of a cstr - + Parameters ---------- c : np.ndarray @@ -114,7 +114,7 @@ def calculate_residual_press_easy_def( """ Calculates the residual equations fo a dead end filtration equation for the pressure in the easy model. - + Parameters ---------- V_dot_P : np.ndarray @@ -138,34 +138,11 @@ def calculate_residual_press_easy_def( return -V_dot_P + deltap * A *hyd_resistance -def calculate_residual_perm_easy_def( - Q_in : float, - V_dot_C : float, - V_dot_P : float - ) -> float: - """ - Calculates the residual equations fo a dead end filtration equation for the permeate Volume - in the easy model. - - Parameters - ---------- - Q_in : float - Flow entering the unit operation - V_dot_P : float - FLow of the Permeate through the membrane and Cake - V_dot_C : float - Volume of the Retentate becoming the Cake - """ - - - return -Q_in + V_dot_C + V_dot_P - def calculate_residual_visc_def(): - """ Calculates the residual of the Viscosity equation of the CSTR """ warnings.warn("Viscosity of def not yet implemented") - return 0 \ No newline at end of file + return 0 diff --git a/CADETPythonSimulator/system_solver.py b/CADETPythonSimulator/system_solver.py index fbd04ca..ce9d6f3 100644 --- a/CADETPythonSimulator/system_solver.py +++ b/CADETPythonSimulator/system_solver.py @@ -423,6 +423,7 @@ def couple_unit_operations( origin_port = origin_info['port'] y_origin_unit = y_initial[self.unit_slices[origin_unit]] + s_unit = origin_unit.get_outlet_state(y_origin_unit, origin_port) s_new += s_unit * Q_destination # Accumulate weighted states diff --git a/CADETPythonSimulator/unit_operation.py b/CADETPythonSimulator/unit_operation.py index a57ce34..6013593 100644 --- a/CADETPythonSimulator/unit_operation.py +++ b/CADETPythonSimulator/unit_operation.py @@ -14,12 +14,11 @@ from CADETPythonSimulator.exception import NotInitializedError, CADETPythonSimError from CADETPythonSimulator.state import State, state_factory from CADETPythonSimulator.residual import ( - calculate_residual_volume_cstr, - calculate_residual_concentration_cstr, + calculate_residual_volume_cstr, + calculate_residual_concentration_cstr, calculate_residual_visc_cstr, calculate_residual_press_easy_def, calculate_residual_cake_vol_def, - calculate_residual_perm_easy_def, calculate_residual_visc_def ) from CADETPythonSimulator.rejection import RejectionBase @@ -115,7 +114,7 @@ def y(self, y: np.ndarray) -> NoReturn: @property def state_derivatives(self) -> dict[str, State]: - """dict: State derivative array blocks of the unit operation, indexed by name.""" + """dict: State derivative array block of the unit operation, indexed by name.""" if self._state_derivatives is None: raise NotInitializedError("Unit operation state is not yet initialized.") @@ -573,12 +572,14 @@ def compute_residual( Q_in = self.Q_in[0] Q_out = self.Q_out[0] - # for i in range(self.n_comp): - # self.residuals['bulk']['c'][i] = c_dot[i] * V + V_dot * c[i] - Q_in * c_in[i] + Q_out * c[i] - # Alternative: Can we vectorize this? - self.residuals['bulk']['c'] = calculate_residual_concentration_cstr(c, c_dot, V, V_dot, Q_in, Q_out, c_in) - self.residuals['bulk']['Volume'] = calculate_residual_volume_cstr(V, V_dot, Q_in, Q_out) + self.residuals['bulk']['c'] = calculate_residual_concentration_cstr( + c, c_dot, V, V_dot, Q_in, Q_out, c_in + ) + + self.residuals['bulk']['Volume'] = calculate_residual_volume_cstr( + V, V_dot, Q_in, Q_out + ) self.residuals['inlet']['viscosity'] = calculate_residual_visc_cstr() @@ -600,34 +601,38 @@ class DeadEndFiltration(UnitOperationBase): """ cake = { 'dimensions': (), - 'entries': {'c': 'n_comp', 'viscosity': 1, 'pressure': 1, 'cakevolume': 1, 'permeate': 1}, + 'entries': {'c': 'n_comp', + 'viscosity': 1, + 'pressure': 1, + 'cakevolume': 1, + 'permeate': 1 + }, 'n_inlet_ports': 1, } - bulk = { + permeate = { 'dimensions': (), 'entries': {'c': 'n_comp', 'viscosity': 1, 'Volume': 1}, 'n_outlet_ports': 1, } - _state_structures = ['cake', 'bulk'] + _state_structures = ['cake', 'permeate'] membrane_area = UnsignedFloat() membrane_resistance = UnsignedFloat() specific_cake_resistance = UnsignedFloat() - molar_volume = SizedUnsignedNdArray(size = 'n_comp') - efficency = SizedUnsignedNdArray(size = 'n_comp') +# molar_volume = SizedUnsignedNdArray(size = 'n_comp') => component system erben von cadet process + rejection = Typed(ty=RejectionBase) _parameters = [ 'membrane_area', 'membrane_resistance', 'specific_cake_resistance', - 'molar_volume', - 'efficency' + 'rejection' ] def compute_residual( self, t: float, - ) -> NoReturn: + ) -> NoReturn: Q_in = self.Q_in[0] Q_out = self.Q_out[0] @@ -643,11 +648,11 @@ def compute_residual( viscosity_in = self.states['cake']['viscosity'] - c = self.states['bulk']['c'] - c_dot = self.state_derivatives['bulk']['c'] + c = self.states['permeate']['c'] + c_dot = self.state_derivatives['permeate']['c'] - V = self.states['bulk']['Volume'] - V_dot = self.state_derivatives['bulk']['Volume'] + V = self.states['permeate']['Volume'] + V_dot = self.state_derivatives['permeate']['Volume'] deltap = self.states['cake']['pressure'] @@ -660,14 +665,53 @@ def compute_residual( # Handle inlet DOFs, which are simply copied to the residual self.residuals['cake']['c'] = c_in - self.residuals['cake']['cakevolume'] = calculate_residual_cake_vol_def(Q_in, efficency, molar_volume, c_in, V_dot_C) - self.residuals['cake']['pressure'] = calculate_residual_press_easy_def(Q_p, V_C, deltap, membrane_area, viscosity_in, membrane_resistance, specific_cake_resistance) - self.residuals['cake']['permeate'] = calculate_residual_perm_easy_def(Q_in, V_dot_C, Q_p) + self.residuals['cake']['cakevolume'] = calculate_residual_cake_vol_def( + Q_in, + efficency, + molar_volume, + c_in, + V_dot_C + ) + + self.residuals['cake']['pressure'] = calculate_residual_press_easy_def( + Q_p, + V_C, + deltap, + membrane_area, + viscosity_in, + membrane_resistance, + specific_cake_resistance + ) + + self.residuals['cake']['permeate'] = calculate_residual_volume_cstr( + V_C, + V_dot_C, + Q_in, + Q_p + ) + self.residuals['cake']['viscosity'] = calculate_residual_visc_def() - self.residuals['bulk']['c'] = calculate_residual_concentration_cstr(c, c_dot, V, V_dot, Q_p, Q_out, c_in) - self.residuals['bulk']['Volume'] = calculate_residual_volume_cstr(V, V_dot, Q_p, Q_out) - self.residuals['bulk']['viscosity'] = calculate_residual_visc_cstr() + new_c_in = (1-efficency)*c_in + + self.residuals['permeate']['c'] = calculate_residual_concentration_cstr( + c, + c_dot, + V, + V_dot, + Q_p, + Q_out, + new_c_in + ) + + self.residuals['permeate']['Volume'] = calculate_residual_volume_cstr( + V, + V_dot, + Q_p, + Q_out + ) + + self.residuals['permeate']['viscosity'] = calculate_residual_visc_cstr() @@ -795,4 +839,4 @@ def compute_residual( Residual of the unit operation. """ - raise NotImplementedError() \ No newline at end of file + raise NotImplementedError() diff --git a/pyproject.toml b/pyproject.toml index 523485d..8aa8299 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -68,4 +68,8 @@ version = { attr = "CADETPythonSimulator.__version__" } [tool.ruff] # Same as Black. line-length = 88 -indent-width = 4 \ No newline at end of file +indent-width = 4 + +[tool.ruff.lint] +select = ["E", "F", "W"] +ignore = ["F401"] \ No newline at end of file diff --git a/tests/test_residual.py b/tests/test_residual.py index 61f3d99..640bd82 100644 --- a/tests/test_residual.py +++ b/tests/test_residual.py @@ -1,13 +1,14 @@ from CADETPythonSimulator.residual import ( - calculate_residual_volume_cstr, calculate_residual_concentration_cstr + calculate_residual_volume_cstr, + calculate_residual_concentration_cstr, + calculate_residual_cake_vol_def, + calculate_residual_press_easy_def ) from CADETPythonSimulator.exception import CADETPythonSimError import pytest import numpy as np -""" -q_in q_out, was sind physikalisch sinnvolle szenarien -""" + # random number test TestCaseConc_level1 = { @@ -37,7 +38,7 @@ "expected" : np.array([0,]) } -# flow in and out are equal, but concentrations going into the unit is not +# flow in and out are equal, but concentrations going into the unit is not TestCaseConc_diffcin = { "values" : { "c" : np.array([0.1,]), @@ -51,7 +52,7 @@ "expected" : np.array([-0.1,]) } -#flow in and out are not equal, concentrantions going in are +#flow in and out are not equal, concentrantions going in are TestCaseConc_diffvol = { "values" : { "c" : np.array([0.1,]), @@ -83,8 +84,8 @@ @pytest.mark.parametrize( "parameters", [ - TestCaseConc_level1, - TestCaseConc_equal, + TestCaseConc_level1, + TestCaseConc_equal, TestCaseConc_diffcin, TestCaseConc_diffvol, TestCaseConc_diffvolc @@ -95,7 +96,10 @@ def test_calculation_concentration_cstr(self, parameters): param_vec_conc = parameters["values"].values() - np.testing.assert_array_almost_equal(calculate_residual_concentration_cstr(*param_vec_conc), parameters["expected"]) + np.testing.assert_array_almost_equal( + calculate_residual_concentration_cstr(*param_vec_conc), + parameters["expected"] + ) @@ -172,7 +176,21 @@ def test_calculation_cstr(self, parameters): param_vec_volume = parameters["values"].values() - np.testing.assert_equal(calculate_residual_volume_cstr(*param_vec_volume), parameters["expected"]) + np.testing.assert_equal( + calculate_residual_volume_cstr(*param_vec_volume), + parameters["expected"] + ) + + + + +class TestResidualCakeDEF(): + def test_calculation_def(self, parameters): + + param_vec_volume = parameters["values"].values() + + np.testing.assert_equal(1,1) + TestCaseConcError = { "values" : { @@ -211,3 +229,4 @@ def test_calculation_concentration_cstr_error(self, parameters): with pytest.raises(CADETPythonSimError): calculate_residual_concentration_cstr(*param_vec_volume) + diff --git a/tests/test_unit_operation.py b/tests/test_unit_operation.py index b4fad0f..fc592d3 100644 --- a/tests/test_unit_operation.py +++ b/tests/test_unit_operation.py @@ -378,9 +378,18 @@ def test_initialize(self, unit_operation: UnitOperationBase, expected: dict): 'Q_out' : [4] }, [ - ("calculate_residual_concentration_cstr", lambda c, c_dot, V, V_dot, Q_in, Q_out, c_in: c_dot * V + V_dot * c - Q_in * c_in + Q_out * c), - ("calculate_residual_visc_cstr", lambda *args : 0), - ("calculate_residual_volume_cstr", lambda V, V_dot, Q_in, Q_out: V_dot - Q_in + Q_out) + ("calculate_residual_concentration_cstr", + lambda c, c_dot, V, V_dot, Q_in, Q_out, c_in: + c_dot * V + V_dot * c - Q_in * c_in + Q_out * c + ), + ("calculate_residual_visc_cstr", + lambda *args : + 0 + ), + ("calculate_residual_volume_cstr", + lambda V, V_dot, Q_in, Q_out: + V_dot - Q_in + Q_out + ) ], { 'inlet' : { @@ -444,9 +453,12 @@ def test_unit_residual( for unit_module, module_dict in expected.items(): for property, value in module_dict.items(): - np.testing.assert_equal(value, unit_operation.residuals[unit_module][property]) + np.testing.assert_equal( + value, + unit_operation.residuals[unit_module][property] + ) # %% Run tests if __name__ == "__main__": - pytest.main(["test_unit_operation.py"]) \ No newline at end of file + pytest.main(["test_unit_operation.py"])