Water quality simulation with bulk coefficient and limiting potential not working properly

[belottif]

Summary
I'm running a water quality simulation using as a test the Net3 network. I added one source with strength 7.4 at the tank '1'. I set a bulk coefficients for each pipe in the network and a global limiting potential of 15. I would expect that the quality at some nodes in the network would increase until reaching 15, but this is not happening. In some cases the quality also decrease instead of increasing. After trying different limiting potentials and different bulk coefficients the results change sometimes based on the bulk coefficients, but it doesn't seem that the limiting potential plays a role. Moreover, for a first order reaction in which unit should I define the bulk coefficients? I guess I should write them in 1/d and then they will transformed by WNTR into 1/s during calculation, is it right?

Example
import wntr
wn = wntr.network.WaterNetworkModel('Net3.inp')
wn.options.quality.parameter = 'CHEMICAL'
wn.add_source('S', '1', 'CONCEN', 7.4,'SourcePattern')
wn.options.reaction.limiting_potential = 15
wn.options.reaction.bulk_coeff = 0.1
lake = wn.get_node('Lake')
lake.initial_quality = 7.4
river = wn.get_node('River')
river.initial_quality = 7.4
for name, junction in wn.junctions():
junction.initial_quality = 7.4

sim = wntr.sim.EpanetSimulator(wn)
base = sim.run_sim()
results = base.node['quality']

Environment
Provide information on your computing environment.

• Operating system: Windows 10 Enterprise
• Python version: 3.11.5
• WNTR version: 1.1.0

[ucchejbb]

@belottif The code snippet below seems to behave the way I would expect in EPANET when I tested it. Please let us know if you have additional questions.

wn = wntr.network.WaterNetworkModel(‘Net3.inp')
wn.options.quality.parameter = 'CHEMICAL'

wn.add_source('S', '1', 'CONCEN', 7.4e-3, 'SourcePattern')
wn.options.reaction.limiting_potential = 15 ## Not currently converted by WNTR
wn.options.reaction.bulk_coeff = 1/(3600*24) ## 1 /day would be input to EPANET, /sec
lake = wn.get_node('Lake').initial_quality = 7.4e-3 ## kg/m**3 or g/L
river = wn.get_node('River').initial_quality = 7.4e-3
for name, junction in wn.junctions():
wn.get_node(junction).initial_quality = 7.4e-3 #7.4
# junction.initial_quality = 7.4 ## this should also work, but I left the alternate way I used

for tank in wn.tank_name_list:
wn.get_node(tank).initial_quality = 7.4e-3

sim = wntr.sim.EpanetSimulator(wn)
base = sim.run_sim()
results = base.node['quality']

for i in ['1', '2', '3', '117']:
plt.plot(results.index/3600, ## to convert to hr
results[i].values* 1e3, ## to convert to mg/L
label=i)

plt.ylabel('Chemical (mg/L)')
plt.xlabel('Time (hours)')
plt.legend()

[belottif]

@ucchejbb Thanks a lot for the help! I can see it was mainly a matter of units, especially the initial quality and the limiting potential still having different conversion.

Just a last question: I would like to set a different limiting potential for each singular pipe. I can see it is possible to code it in WNTR, but it is kind of ignored during calculation. This is probably due to EPANET, but as far as I know it is possible to do that with EPANET-MSX. Do you plan to implement this option also in WNTR?

Thanks for your reply!

[dbhart]

Hi @belottif
EPANET-MSX extension is currently in pull request #426
and this will be pulled into WNTR 1.3. I'm hopeful that will happen very soon now that tests are passing (this also fixes many of the macos ARM vs x86 issues). If you want to follow that PR you will get notified when it gets closed and pulled into main.