MAW head not correctly modeled when its head is below screen bottom

[ougx]

When the head of GWF nodes connected to the well falls below the well bottom, the model wont update the MAW head correctly.

If the head of the connected cell is below a MAW well bottom, the MAW head wont be properly updated to represent dry condition.

An single cell and single well example showing this problem can be downloaded at: https://github.com/MODFLOW-USGS/modflow6/files/5068947/MAW_Unit_Test.zip

[ougx]

I am new to the mf6 code and trying to tracking down this problem. I am trying to find the piece of source code that updates the MAW head after the solution.

Is it the MAW head associated by NumericalSolutionType%x and will be updated automatically when the Ax=b is solved? Or are there some lines updating the MAW head such as mawType%head = NumericalSolutionType%x(i1:i2)?

Thanks!

[langevin-usgs]

When you say that the "the model wont update the MAW head correctly" what MAW head did you get and what MAW head did you expect?

MawType%xnewpak is a pointer and points into the correct range within NumericalSolutionType%x. So MawType%xnewpak will always be the latest solution from the solver for the multi-aquifer wells in the package.

I tried to download your test model but the link is not active. Please zip up the model input files and drag them into this issue, and we'll take a look.

[ougx]

MAW_Unit_Test.zip

See attached. The well bottom is at 15 but the head of the groundwater node is 4.5. When NumericalSolutionType%x(i_maw) was updated to 4.5 in the solver, somehow it will increase to some other values (and I don't see the pattern of these new values) in maw_calculate_conn_terms. Below is an example output (the fourth node is the MAW). The head of MAW is calculated to be 4.5 as the head of the connected GW nodes. However, it becomes 12.xx in maw_calculate_conn_terms (no backtracking or under relaxation).

Outer Inner:           9           2
                                                   solve(): this@xtemp  4.5000000      4.5000000      4.5000000      31.531390    
                                                       solve(): this@x  4.5000000      4.5000000      4.5000000      4.5000000    
                     maw_calculate_conn_terms: this@xnew, this@xnewpak  4.5000000      4.5000000      4.5000000      12.153139      12.153139   

[ougx]

Looks like the MAW head is updated in maw_nur. It used the well bottom for zmin. It will bring the well alive every time the well head is drop below the well bottom, though the connected cells are still dry. I believe this value should be the lowest cell's bottom elevation, which is smaller than well bottom in this example.

[langevin-usgs]

Hi @ougx, just wanted to let you know that we are looking into this. I think we are slowly getting to the bottom of it and hope to have a fix within the next few days.

[ougx]

Thanks. Another issue that may need attention is that a large pumping rate which will make the well dry will create convergence issue. See the attached example with the similar setup. And when I tried to set HEAD_LIMIT at different heights, I got different outcomes:

# 1 head_limit 10 ### floating overflow 
# 1 head_limit 11 ### runs but not converge
# 1 head_limit 12 ### runs and converge

MAW_Unit_Test2.zip

[langevin-usgs]

See #549 for a MAW fix that should handle the problems you've been having.

[langevin-usgs]

@ougx we are finishing up this bug fix with #549. Regarding your second issue with the HEAD_LIMIT option, we are not entirely sure this is a problem. As indicated in the mf6io.pdf document, the RATE_SCALING option is a better way to smoothly reduce the pumping rate in a well as the well head approaches a limit. If you continue to look into the HEAD_LIMIT option and find a bug in the code, feel free to submit a pull request. Thanks again for your help identifying this issue with MAW.