use partitioning on real-world model

[luitjansl]

The partitioning of flow and transport simulations has been tested on simple models, but we would like to try it on a more complex, real-world simulation. For that we can use Teun's model, which can be found here: n:\Deltabox\Postbox\Engelen,van Joeri\test_mf6\

[JoerivanEngelen]

Update

Splitting the model worked succesfully, after fixing two validation mistakes. However, comparing heads shows differences when the buoyancy package is active. I haven't found the precise cause of this out.

The following plots show the differences:

(head["head"] - head_split["head"]).isel(time=ind[0], layer=10)

With buoyancy package:

Without buoyancy package:

With buoyancy package, but density slope = 0.0, so that density has no effects on heads but the reformulation is active:

Things I have tried that had no effect before investigating buoyancy package.

• In the exchange file, change IHC=1 to IHC=2. This even worsened errors with a factor 4.
• Remove vertical staggering of grid, by setting all top and bottoms to one single value (instead of differing over x,y). No effect.
• Set the solver settings of both the transport as well as flow solver stricter (reducing DVCLOSE by a factor 10). No effect.
• Remove advection package. No effect.

[JoerivanEngelen]

Update 2:

I tried loading the model with Flopy and split it, which took some work as Flopy could not read the griddata as floats (I think because of the NAN's we write), so I had to simplify all griddata with NaNs to layered constants, write as absolute paths, and write textfiles instead of binary. Consequently, I could split the model, but the latest released Flopy version (3.6) didn't support splitting transport models yet. However, the gwfgwf files it produced appeared to contain the same values, also auxiliary, so I don't think this is going to lead us somewhere else.

Furthermore, I reduced the amount of topsystem packages to just a GHB and RCH package, differences changed, but didn't reduce really.

[JoerivanEngelen]

Update 3:

Contacted @mjr-deltares to see if this might be a Modflow 6 issue, as multi-model variable-density flow models have not been tested extensively (I think). Closing this for now, will re-open if it turns out this is an iMOD Python issue.