PEST-MOU - Pumping Optimization

[gmerritt123]

I'm trying to set up an optimization of maximizing pumping rates without letting water levels in wells dip draw down belows certain levels. Quick rundown:

1. decision variables: assigned pumping rates at ~45 wells

2. uncertain/adjustable model parameters: none (for now)

3. The model itself is rough surrogate "response" model kinda based on https://www.esaa.org/wp-content/uploads/2015/01/WaterTech2014-p48.pdf . It takes model runtime from ~10 mins to ~1 second and makes (approximate) optimization computationally feasible

4. Observations coming out of the model:

• "Qw", just the assigned rate at each well. These are listed in mou_objectives, all with weight of 1, wanting to max this
• "g_dd", simulated remaining drawdown, set as a greater than constraint, e.g. "well water level can't go below 10 cm"
• "wf_cons", certain combinations of pumping wells can only pump at a certain combined maximum due to treatment capacity. E.g. "wells 1-3 combined pumping can't exceed 120 L/s."

One kinda big wrinkle, similar to knapsack problem: if a well isn't pumping, we don't care if it goes dry. I initially formulated this in the post processing of the model: "if applied pumping rate is 0, set its drawdown to exactly 10 cm (i.e. guaranteed constraint met). I realized this won't really work, because since the parameter (pumping) bounds are continuous between 0 and x, a population member actually pumping exactly 0 at it is nigh impossible to happen, and even if does, the odds that it would combine with a better "optimal" solution that pest would see is even less likely. So what I've done for now is instead set a threshold, saying "if pumping is less than 3 L/s, set constraint = met". This is also imperfect...

I'm considering implementing (discussed in another thread) an "on-off" switch decision for each well, set to between 0 and one, and if <0.5 then hard set rate to 0 and meet constraint". This kinda messes with my gradient-based solver brain though, because pest could assign a sky high pumping rate but have the switch off, and that kinda gives it conflicting information no? Is there a parameter tying trick to this? Any thoughts/advice on this appreciated...

For initial populations, I've tried both uniform (default) initial population and a "tighter" gaussian draw around an initial "eyeball optimized" guess. I've used DE and PSO generators.

What is happening on my current runs:

1. Over the first few generations, it seems to be really really good at getting feasible solutions (no constraint violations), but...
2. After literally hundreds of generations, it's still staying feasible but not finding anything better than my initial guess! With every passing generation it seems to be "stuck" with total pumping hanging around 90% of my initial best guess.

Anything diagnostic I can look at? Because it's a surrogate model, the whole model package is literally a python script (reliant on numpy and pandas,) two pickled dataframes, and 1-2 lookup csvs. So the entire thing is <80 kb (~ 4 mb if you count the pest executeable). So I've attached it all zipped up here if anyone feels like taking a gander. As always, thanks.

v1_mou.zip

[cnicol-gwlogic]

A couple of things:

1. The on/off thing works, I've done it. Similar optimisation problem to yours.
2. My understanding is that the algorithms aren't so much gradient based, they're more selection based (not that I've read up on the details...). Regardless, it worked for me.
3. Re stymied progress - is it just possible that you've got your answer? You mention quite a few generations, so I wouldn't be surprised if it got there. You could do some manual checks - tweak a well rate or two from a selected member (up), and see what happens to your 10cm head constraint. And/or just manually run a selected optimal member and inspect the heads and rates etc - is the workflow spitting out the correct level / flow rate info for mou to read? Does the manual run make sense?
4. With the on/off option, I have also included a less than objective on number of operating wells (cost / infrastructure minimisation).
5. One thing I have had to be wary of with these types of MoU runs is the wel package i use. I usually use cln / maw style wells with defined well flow cutoff levels (minima), because well bore storage / losses and pump / screen depths etc are important for borefield reliability (max yield / min level) type work, particularly when looking at extremes like drought yields. If using these types of wells (with level based flow cutoffs), rather than having an obs of applied flow rate, you need an obs of achieved flow rate (eg from the Cbb file). May or may not be relevant for you.
6. If concerned about the stymied progress still, you could try sample your initial well rates down the bottom end of the range, then let mou push things up gradually to map out the front.

[jtwhite79]

In addition to @cnicol-gwlogic 's good advice (which is always good!)...I think there might be a conceptual issue with that many objectives. There is a distinction between multiobjective and "many objective", where the former is well understood and is usually a max of 4-5 objectives. Many objective opt is still an open problem AFAIK. While the nsga-II alg as coded in mou with mechanically handle unlimited objectives, conceptually when you have heaps of objectives, almost every member ends up being "non-dominated" (in the nsga-II definition) so that the alg cant really evolve preferred characteristics - so in essence, all the beta males are reproducing haha!

Any chance you can aggregrate those objectives into 1 or 2 objectives? You dont have to worry about nonlinearity/thresholds/discontinuities in the aggregated values...

[cnicol-gwlogic]

That's true, sorry I hadn't looked at the script. I've found even like 5 objectives gets hairy.

[gmerritt123]

Thanks both of you for advice and taking a look.

Yes, the way I formulated the objectives was also a big Q I didn't mention, I was not confident that my "max each well rate, weight them all the same" formulation would give the right info to the evolution algos.... so yeah, I can literally distill it down to one objective (maximize the sum of all pumping). I did try this in one iteration/setup I tried... and it didn't seem to help... but there could have been other things in that particular iteration and now that you've chimed in I now know to go back to it. So thanks again for taking a look.

I think my next setup based on your advice will be 1) go to single objective, max sum of all pumping, and 2) implement the on-off switch, and 3) play with initial conditions. The real uncertainty however is when I'll have time to revisit this thing as it's firmly in side-project territory at this point. I'll report back here whenever I do. Thanks again!