Pestpp-da advice

[laat0003]

Hi,

I'm working on a project where I have a fair amount of transient data for relatively simple MF6 model (single layer unconfined, quadtree grid, olproc, plproc with pilot points, seglists and gplanes). I already have it set up to run with IES but am keen to try PESTPP-DA specifically because the area is influenced by multiple mining operations so initial conditions for the model are unknown. I have two questions. If I am estimating the initial conditions for each ensemble member, then a steady-state stress period at the start of the model is not required, correct? The other question is to do with cycle lengths. Are they completely arbitrary, or is there a recommended strategy? Here I am thinking more along the lines of rules of thumb like "more is better than less", "break periods with high observation density into separate cycles" that sort of thing.

Thanks in advance for your time.
Tariq

[jtwhite79]

Hey @laat0003 - sounds like an interesting use case for pestpp-da. Conceptually, the choice around how to setup cycles is multifaceted - it needs to balance the frequency of the observations, the frequency/temporal character of the predictions (short-term operational vs long-term), as well as computational burden (so in some way, the same as it ever was!). I wonder if you could estimate initial conditions with ies if you are worried about their impact? If you are only worried about getting the initial conditions "on track" maybe just a few short-term cycles followed by one long cycle would help with this? In the sequential framework, I dont think you want an initial SS period because that would potentially "bake in" SS conditions that then the subsequent cycles would have to work to overcome/assimilate (unless those SS conditions are an acceptable "first order approximation" and then you want to use the next few cycles to transition from SS to transient)

If you do go forward with sequential joint state-parameter estimation, a few things to consider:

• you'll need to add the final simulated water level at every node as an "observation" and you'll need to add the initial water level for every node as a parameter - these two quantities will then need to be "linked" either thru identical naming or by adding a special column to the observation data
• Ive found its likely easier to sample the observed data to the cycle ends and use the obs cycle table and weight cycle table for observations that are also the states you are estimating (in this case, the groundwater levels), so it might be easier to pull those obs out of olproc and just get them from the raw sim water levels - I realize this ignores spatial and/or temporal interpolation but the alternative might include needing a dedicated olproc run for each cycle...
• the computational burden of sequential estimation is higher that the smoother so just be ready for that!

@aymanalz will likely have some more insightful things to say...

[aymanalz]

I agree with you , @jtwhite79. In principle, data assimilation can update the initial head (or any other simulated state) and produces statistically sound initial head field; however, the resulting initial head might not be physically correct. The difference might be small, especially when observations are clean of errors and cover most of simulated domain, but we should keep this in mind. If transient runs have a steady-state run at the beginning, then the initial head is guaranteed to honor the physics of groundwater flow.

In my opinion, updating the head sequentially can be very useful when the model itself has structural errors (e.g. missing pumping or incorrect boundary condition, or any model issue that produces biased predictions). In this case, PESTPP-DA will try to eliminate this bias every cycle and restarts the model using 'correct' initial head. It's very important to notice here that when the model has structural errors, then PESTPP-DA will try compensate for that by producing extreme parameters values (e.g. very high (or small) conductivity values). PEST++ family (including PESTPP-DA) controls this situation by imposing upper and lower bounds on parameter values. When a parameter reaches one of its bounds and model predictions are still has large deviation from observations (underfitting), then PESTPP-DA might update model predictions (e.g. simulated head) to correct for errors that weren't removed by changing the unknown parameters. Sorry for my long comment.

[laat0003]

Thanks to you both @aymanalz and @jtwhite79 for your very informative responses, which have provided me with some much needed insight. I had considered the short initial cycle and then one long one as a suitable strategy in order to get the initial heads on track but needed some confirmation that this was indeed worthwhile. Truth be told, it was trying to reconcile the use of olproc that prompted my questions to begin with so I am grateful for the suggestions from @jtwhite79. No need for any apologies regarding length of responses @aymanalz, it is all very useful. Thanks again.