Updating EN_LEAK_AREA during simulation produces inconsistent EN_LINK_LEAKAGE and EN_LEAKAGEFLOW behavior

[dapaez270]

When updating EN_LEAK_AREA during an EPS hydraulic simulation using ENsetlinkvalue(), the leak area value appears to update successfully (as returned by ENgetlinkvalue()), and the reported EN_LINK_LEAKAGE / EN_LEAKAGEFLOW values also change accordingly. However, the hydraulic solution itself does not seem to reflect the updated leak condition.

Specifically, after setting EN_LEAK_AREA = 0 during the simulation, the reported leakage becomes zero, but the total pipe flow (EN_FLOW) remains unchanged, suggesting that the network hydraulics are still being solved using the previous leak configuration.

Example:

See attached model and following script (python wrapper based on owa-epanet):

inp_file = r"One branch.inp"
rpt_file = inp_file[:-4] + ".rpt"
out_file = inp_file[:-4] + ".out"

em.ENopen(inp_file, rpt_file, out_file)

id_pipe = em.ENgetlinkindex("C2")

em.ENsetlinkvalue(id_pipe, em.EN_LEAK_AREA, 30)

em.ENopenH()
em.ENinitH(flag=em.EN_SAVE)
t = 0
print(f"Time(s) \t LeakA(mm2) \t LeakFlow(Lps) \t TotFlow (Lps)")

while True:
    
    if t < 6*60*60: #6hrs into simulation
        em.ENsetlinkvalue(id_pipe, em.EN_LEAK_AREA, 30)
    else:
        em.ENsetlinkvalue(id_pipe, em.EN_LEAK_AREA, 0)
    
    t = em.ENrunH()
    
    leak_area = em.ENgetlinkvalue(id_pipe, em.EN_LEAK_AREA)
    leak_flow = em.ENgetlinkvalue(id_pipe, em.EN_LINK_LEAKAGE)
    tot_flow = em.ENgetlinkvalue(id_pipe, em.EN_FLOW)
    
    print(f"{t} \t {leak_area} \t\t {leak_flow:.4f} \t\t {tot_flow:.4f}")
    
    tstep = em.ENnextH()
    if tstep <= 0:
        break
        
em.ENcloseH()
em.ENsaveH()   
em.ENclose()  

Output

Time(s) 	 LeakA(mm2) 	 LeakFlow(Lps) 	 TotFlow (Lps)
0  		 30.0 		 1.5828 		 5.7869
3600 	 30.0 		 1.5828 		 5.7869
7200 	 30.0 		 1.5828 		 5.7869
10800 	 30.0 		 1.5828 		 5.7869
14400 	 30.0 		 1.5828 		 5.7869
18000 	 30.0 		 1.5828 		 5.7869
21600 	 30.0 		 1.5828 		 5.7869
25200 	 0.0 		 0.0000 		 5.7869
28800 	 0.0 		 0.0000 		 5.7869
32400 	 0.0 		 0.0000 		 5.7869
...
86400 	 0.0 		 0.0000 		 5.7869

Note that at t = 25200 the EN_LEAK_AREA is zero, and the reported EN_LINK_LEAKAGE is zero, but the EN_FLOW remains as 5.7869 meaning there is still a computed leakage.

I think the behaviour of EN_LEAK_AREA should be the same as other properties like EN_DIAMETER which support being updated during the EPS simulation, but if not, the ENgetlinkvalue() and ENgetnodevalue() functions should return the actual unchanged values for EN_LEAK_AREA, EN_LINK_LEAKAGE, and EN_LEAKAGEFLOW.

One branch.inp.txt

[LRossman]

Before a simulation begins, a pipe's leak area parameter (Ao) is converted into an emitter coefficient for both of its end nodes which is added to the coefficients from other leaking pipes connected to the nodes:

Ca = 1 / ($\sum$ (4.815 Ao))^2

where Ca is the node's overall emitter coefficient for leak area and the summation is over all leaking pipes connected to it. The actual leakage rate as a function of pressure is computed using the nodal emitter. Once hydraulic convergence is achieved, the leakage rate for the pipe is computed directly using the pipe's Ao value.

If a pipe's Ao changes at some point during a simulation then Ca would have to be updated to reflect the change and influence nodal leakage. This currently isn't done in the code. A similar updating to the emitter coefficient for leak expansion would have to be made if the leak expansion parameter were changed.

[dapaez270]

Thank you, Lew. This clarifies the current behavior.

It sounds like dynamic updating of leak area and leak expansion during a simulation would fall more under a potential enhancement for future versions. At the same time, the reporting of leakage flows based on attempted updates to these variables may be considered a minor inconsistency/bug, since the reported values appear to reflect changes that are not actually being applied hydraulically.

Probably not something that needs immediate attention, as others trying to do this may find this thread, but hopefully something to consider for future releases.

[LRossman]

A simple fix for the following code section in EN_setlinkvalue:

    case EN_LEAK_AREA:  // leak area per 100 pipe lengths units
        if (value < 0.0) return 211;
        Link[index].LeakArea = value * Ucf[LENGTH];
        if (Hyd->OpenHflag) openleakage(p);   //<------ new line
        break;

with a similar fix for EN_LEAK_EXPAN. The openleakage() function is an existing function in leakage.c that is called at the start of a hydraulic analysis to allocate and initialize leakage coefficients for the end nodes of pipes that leak.

This is a somewhat inefficient solution since it frees and then re-allocates and initializes all of the node leakage coefficients if just a single pipe has its leakage parameter changed. I tried the more efficient approach of adjusting the current set of node coefficients for just the pipe's end nodes but it took over 100 lines of code to do so and I couldn't get it to work for the case where a pipe leakage parameter was changed to 0.

If no one objects I'll submit a PR for this solution.

[LRossman]

Here's a revised fix for setting leakage parameters while a hydraulic simulation is underway that is more efficient than the previous one:

1. Add a LeakageChanged flag to the SLeakage struct in types.h.
2. In the EN_setlinkvalue function of epanet.c, set the flag to TRUE when a pipe's leakage area or expansion is changed while the OpenHflag is TRUE (meaning a hydraulic analysis is underway).
3. In the runhyd function of hydraul.c, call the openleakage function to reset the leakage system if LeakageChanged is TRUE.
4. In the openleakage function of leakage.c, re-set LeakageChanged to FALSE to prevent unneccessary calls of it in subsequent time periods.

I have confirmed that these changes, included in a new dev-set_leakage_fix branch, work for the single leaking pipe example of @dapaez270 . I hope to do additional testing on more complex networks and invite others to also do so. If testing proves successful I will submit a PR to merge the branch into dev.

[dapaez270]

I tested the dev-set_leakage_fix branch using the single leaking pipe model and a second complex network model with multiple leaking pipes with varying leakage parameters.

The proposed fix appears to resolve the issue on my end.

Thank you for implementing and sharing the branch for testing!