Skip to content

Commit

Permalink
Add molar flows for precipitate and reagent to simplify reaktoro inte…
Browse files Browse the repository at this point in the history 
…gration (#1498)

* Update stoichiometric_reactor.py

* Update stoichiometric_reactor.rst

* fixed tests

* fixed pylint

* update tnameing for consistency

* Apply suggestions from code review

---------

Co-authored-by: Adam Atia <aatia@keylogic.com>
  • Loading branch information
@avdudchenko @adam-a-a
avdudchenko and adam-a-a authored Dec 6, 2024
1 parent d886bd4 commit 8a0c9ed
Show file tree
Hide file tree
Showing 3 changed files with 100 additions and 39 deletions.
2 changes: 2 additions & 0 deletions docs/technical_reference/unit_models/stoichiometric_reactor.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -133,9 +133,11 @@ Variables
"Reagent dose", reagent_dose,[reagent],kg/:math:`\text{m}^3`
"Reagent density", density_reagent,[reagent],kg/:math:`\text{m}^3`
"Reagent flow mass", flow_mass_reagent,[reagent],kg/s
"Reagent molar flow", flow_mol_reagent,[reagent],mol/s
"Reagent flow volume", flow_vol_reagent,[reagent],:math:`\text{m}^3`/s
"Stoichiometric coefficients for dissolution", dissolution_stoich_comp, "[reagent, :math:`j`]",dimensionless
"Flow mass of precipitant",flow_mass_precipitate,[precipitant],kg/s
"Molar flow of precipitant",flow_mol_precipitate,[precipitant],mol/s
"Mass concentration of precipitant",conc_mass_precipitate,[precipitant],kg/:math:`\text{m}^3`
"Stoichiometric coefficients for precipitation", precipitation_stoich_comp, "[precipitant, :math:`j`]",dimensionless
"Fraction of solids in waste stream", waste_mass_frac_precipitate, None, fraction
Expand Down
129 changes: 91 additions & 38 deletions watertap/unit_models/stoichiometric_reactor.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,6 @@
Var,
Param,
Suffix,
NonNegativeReals,
Reals,
Reference,
TransformationFactory,
Expand Down Expand Up @@ -281,42 +280,38 @@ def build(self):
)
# Add unit parameters
if self.has_dissolution_reaction:
self.mw_reagent = Var(
self.mw_reagent = Param(
self.reagent_list,
units=pyunits.kg / pyunits.mol,
doc="Molecular weight of reagents",
mutable=True,
)
self.density_reagent = Var(
self.density_reagent = Param(
self.reagent_list,
initialize=1000,
units=units_meta("mass") / units_meta("volume"),
doc="Density of reagents",
mutable=True,
)
for r in self.reagent_list:
self.mw_reagent[r].fix(
pyunits.convert(
self.config.reagent[r]["mw"],
to_units=pyunits.kg / pyunits.mol,
)
self.mw_reagent[r] = pyunits.convert(
self.config.reagent[r]["mw"],
to_units=pyunits.kg / pyunits.mol,
)
if self.config.reagent[r].get("density_reagent") is not None:
self.density_reagent[r].fix(
pyunits.convert(
self.config.reagent[r].get("density_reagent"),
to_units=(units_meta("mass") / units_meta("volume")),
)
self.density_reagent[r] = pyunits.convert(
self.config.reagent[r].get("density_reagent"),
to_units=(units_meta("mass") / units_meta("volume")),
)
else:
_log.warning(
"User did not specify density for reagent {}, using 1kg/L".format(
r
)
)
self.density_reagent[r].fix(
pyunits.convert(
1 * pyunits.kg / pyunits.liter,
to_units=(units_meta("mass") / units_meta("volume")),
)
self.density_reagent[r] = pyunits.convert(
1 * pyunits.kg / pyunits.liter,
to_units=(units_meta("mass") / units_meta("volume")),
)
self.dissolution_stoich_comp = Param(
self.reagent_list,
Expand All @@ -334,21 +329,28 @@ def build(self):
self.reagent_dose = Var(
self.reagent_list,
initialize=1,
domain=NonNegativeReals,
domain=Reals,
units=units_meta("mass") / units_meta("volume"),
doc="reagent dose",
)
self.flow_mass_reagent = Var(
self.reagent_list,
initialize=1e-3,
domain=NonNegativeReals,
domain=Reals,
units=units_meta("mass") / units_meta("time"),
doc="Mass flowrate of reagent",
)
self.flow_mol_reagent = Var(
self.reagent_list,
initialize=1e-3,
domain=Reals,
units=pyunits.mol / units_meta("time"),
doc="Mass flowrate of reagent",
)
self.flow_vol_reagent = Var(
self.reagent_list,
initialize=1e-3,
domain=NonNegativeReals,
domain=Reals,
units=units_meta("volume") / units_meta("time"),
doc="Volume flowrate of reagent",
)
Expand All @@ -368,17 +370,16 @@ def build(self):
units=pyunits.dimensionless,
doc="Solid mass fraction in sludge",
)
self.mw_precipitate = Var(
self.mw_precipitate = Param(
self.precipitate_list,
units=pyunits.kg / pyunits.mol,
doc="Molecular weight of precipitate",
mutable=True,
)
for p in self.precipitate_list:
self.mw_precipitate[p].fix(
pyunits.convert(
self.config.precipitate[p]["mw"],
to_units=pyunits.kg / pyunits.mol,
)
self.mw_precipitate[p] = pyunits.convert(
self.config.precipitate[p]["mw"],
to_units=pyunits.kg / pyunits.mol,
)

self.precipitation_stoich_comp = Param(
Expand All @@ -399,14 +400,21 @@ def build(self):
self.flow_mass_precipitate = Var(
self.precipitate_list,
initialize=1e-3,
domain=NonNegativeReals,
domain=Reals,
units=units_meta("mass") / units_meta("time"),
doc="Mass flowrate of precipitate",
)
self.flow_mol_precipitate = Var(
self.precipitate_list,
initialize=1e-3,
domain=Reals,
units=pyunits.mol / units_meta("time"),
doc="Mass flowrate of precipitate",
)
self.conc_mass_precipitate = Var(
self.precipitate_list,
initialize=1,
domain=NonNegativeReals,
domain=Reals,
units=units_meta("mass") / units_meta("volume"),
doc="Mass concentration of precipitate",
)
Expand Down Expand Up @@ -585,6 +593,14 @@ def eq_dissolution_reaction_generation(b, t, j):
for r in self.reagent_list
)

@self.Constraint(
self.flowsheet().config.time,
self.reagent_list,
doc="Mole generation from dissolution",
)
def eq_flow_mol_reagent(b, t, r):
return b.flow_mass_reagent[r] / b.mw_reagent[r] == b.flow_mol_reagent[r]

@self.Constraint(
self.flowsheet().config.time,
self.reagent_list,
Expand All @@ -604,8 +620,8 @@ def eq_flow_mass_reagent(b, t, r):
)
def eq_flow_vol_reagent(b, t, r):
return (
b.flow_mass_reagent[r]
== b.flow_vol_reagent[r] * b.density_reagent[r]
b.flow_mass_reagent[r] / b.density_reagent[r]
== b.flow_vol_reagent[r]
)

if self.has_precipitation_reaction:
Expand Down Expand Up @@ -649,6 +665,17 @@ def eq_precipitation_reaction_generation(b, t, j):
for p in self.precipitate_list
)

@self.Constraint(
self.flowsheet().config.time,
self.precipitate_list,
doc="Mole generation from precipitation",
)
def eq_flow_mol_precipitate(b, t, p):
return (
b.flow_mass_precipitate[p] / b.mw_precipitate[p]
== b.flow_mol_precipitate[p]
)

@self.Constraint(
self.flowsheet().config.time,
self.precipitate_list,
Expand Down Expand Up @@ -784,7 +811,10 @@ def calculate_scaling_factors(self):
# these variables should have user input, if not there will be a warning

if self.has_dissolution_reaction:

for r in self.reagent_list:
sf = 1 / self.mw_reagent[r].value
iscale.set_scaling_factor(self.mw_reagent[r], sf)
if iscale.get_scaling_factor(self.reagent_dose[r]) is None:
if iscale.get_scaling_factor(self.flow_mass_reagent[r]) is not None:
# scale reagent_dose based on flow_mass_reagent
Expand All @@ -802,15 +832,10 @@ def calculate_scaling_factors(self):
self.reagent_dose[r], default=1, warning=True
)
iscale.set_scaling_factor(self.reagent_dose[r], sf)
if iscale.get_scaling_factor(self.flow_mass_reagent[r]) is None:
if iscale.get_scaling_factor(self.flow_vol_reagent[r]) is None:
sf = iscale.get_scaling_factor(
self.dissolution_reactor.properties_in[0].flow_vol_phase["Liq"]
self.flow_mass_reagent[r], default=1e4, warning=True
)
sf = iscale.get_scaling_factor(self.reagent_dose[r]) * sf
iscale.set_scaling_factor(self.flow_mass_reagent[r], sf)
if iscale.get_scaling_factor(self.flow_vol_reagent[r]) is None:
sf = iscale.get_scaling_factor(self.flow_mass_reagent[r])
print(sf, self.density_reagent[r].value)
sf = sf / self.density_reagent[r].value
iscale.set_scaling_factor(self.flow_vol_reagent[r], sf)
for (t, j), v in self.dissolution_reaction_generation_comp.items():
Expand Down Expand Up @@ -855,9 +880,23 @@ def calculate_scaling_factors(self):
self.dissolution_reactor.properties_in[0].flow_vol_phase["Liq"]
)
iscale.set_scaling_factor(self.flow_mass_reagent[r], sf)
if iscale.get_scaling_factor(self.flow_mol_reagent[r]) is None:
sf = iscale.get_scaling_factor(self.flow_mass_reagent[r])
sf = sf / iscale.get_scaling_factor(self.mw_reagent[r])
iscale.set_scaling_factor(self.flow_mol_reagent[r], sf)
for (t, r), con in self.eq_flow_mol_reagent.items():
sf = iscale.get_scaling_factor(
self.flow_mol_reagent[r], default=1, warning=True
)
iscale.constraint_scaling_transform(con, sf)
iscale.constraint_scaling_transform(
self.dissolution_reactor.eq_isothermal_dissolution[0], 1e-2
)

if self.has_precipitation_reaction:
for p in self.precipitate_list:
sf = 1 / self.mw_precipitate[p].value
iscale.set_scaling_factor(self.mw_precipitate[p], sf)
if iscale.get_scaling_factor(self.flow_mass_precipitate[p]) is None:
if (
iscale.get_scaling_factor(self.conc_mass_precipitate[p])
Expand All @@ -878,7 +917,21 @@ def calculate_scaling_factors(self):
self.flow_mass_precipitate[p], default=1e3, warning=True
)
iscale.set_scaling_factor(self.flow_mass_precipitate[p], sf)
if iscale.get_scaling_factor(self.flow_mass_precipitate[p]) is None:
iscale.set_scaling_factor(self.flow_mass_precipitate[p], sf)
if iscale.get_scaling_factor(self.flow_mol_precipitate[p]) is None:
sf = iscale.get_scaling_factor(self.flow_mass_precipitate[p])
sf = sf / iscale.get_scaling_factor(self.mw_precipitate[p])

iscale.set_scaling_factor(self.flow_mol_precipitate[p], sf)
iscale.constraint_scaling_transform(
self.precipitation_reactor.eq_isothermal_precipitation[0], 1e-2
)
for (t, p), con in self.eq_flow_mol_precipitate.items():
sf = iscale.get_scaling_factor(
self.flow_mol_precipitate[p], default=1, warning=True
)
iscale.constraint_scaling_transform(con, sf)
for p in self.precipitate_list:
if iscale.get_scaling_factor(self.conc_mass_precipitate[p]) is None:
sf = iscale.get_scaling_factor(
Expand Down
8 changes: 7 additions & 1 deletion watertap/unit_models/tests/test_stoichiometric_reactor.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -462,6 +462,10 @@ def test_molar(self, basic_unit_molar):
m.fs.unit.flow_mass_reagent["CaO"] / (56.0774 * 1e-3)
)

assert pytest.approx(flow_ca_in_reagent, rel=1e-5) == value(
m.fs.unit.flow_mol_reagent["CaO"]
)

assert pytest.approx(expected_ca_mol_flow, rel=1e-5) == value(
m.fs.unit.dissolution_reactor.properties_out[0].flow_mol_phase_comp[
"Liq", "Ca_2+"
Expand Down Expand Up @@ -617,7 +621,9 @@ def test_precipitation(self, precipitation_reactor):
"Liq", "Ca_2+"
]
)

assert pytest.approx(1e-3 / (100.09 * 1e-3), rel=1e-5) == value(
m.fs.unit.flow_mol_precipitate["Calcite"]
)
expected_mg_mol_flow = 0.1 - 1e-4 / (58.3197 * 1e-3)

assert pytest.approx(expected_mg_mol_flow, rel=1e-5) == value(
Expand Down

0 comments on commit 8a0c9ed

Please sign in to comment.