gwmod / nlmod Goto Github PK
View Code? Open in Web Editor NEWPython package to build, run and visualize MODFLOW 6 groundwater models in the Netherlands.
Home Page: https://nlmod.readthedocs.io
License: MIT License
Python package to build, run and visualize MODFLOW 6 groundwater models in the Netherlands.
Home Page: https://nlmod.readthedocs.io
License: MIT License
When running the nlmod.mdims.mlayers.combine_layers_ds
function in the example notebook 04_modifying_layermodels I get this error:
TypeError: Could not convert tuple of form (dims, data[, attrs, encoding]): (9, 10, 11, 12, 13) to Variable.
I think it has something to do with commit 17ed300 because the error is exactly in the modified line.
if fill_nans=True
in the function nlmod.util.get_heads_array
it automatically assumes that the nan_value is the maximum value. Sometimes the nan value is -1.e+30. Make the nan_value an optional argument of this function.
The figure of the recharge in example notebook 3 about local grid refinement show that the recharge is sometimes zero.
I think we should choose a specific grid-type, so we do not have to maintain a lot of code to support multiple grid-types. For almost every calculation some sort of grid-refinement will be needed, for example a coarser grid at the boundaries of the model.
Therefore I would propose to only support the disv-package (so what we call unstructured), even if the grid has constant column width and row height. In this way we can optimize our code for this special case, and generate methods to quickly transform and display irregular data on a regular grid, and vice versa.
When I run the Bergen model twice, not all cached data is reused, as can be seen below.
Intersecting oppwater with grid: 100%|██████████| 9/9 [00:02<00:00, 4.17it/s]
Loading gridded surface water data from cache.
INFO:nlmod.mdims.mlayers:get active cells (idomain) from bottom DataArray
INFO:nlmod.mdims.mgrid:get first active modellayer for each cell in idomain
INFO:nlmod.mdims.mlayers:using top and bottom from model layers dataset for modflow model
INFO:nlmod.mdims.mlayers:replace nan values for inactive layers with dummy value
INFO:nlmod.mdims.mlayers:add kh and kv from model layer dataset to modflow model
INFO:nlmod.mdims.mlayers:nan values at the northsea are filled using the bathymetry from jarkus
INFO:nlmod.cache:cache was created using different numpy array values, do not use cached data
INFO:nlmod.cache:cache was created using different dictionaries, do not use cached data
INFO:nlmod.mfpackages.mfpackages:creating modflow SIM, TDIS, GWF and IMS
INFO:nlmod.cache:caching data -> sea_model_ds.nc
INFO:nlmod.cache:cache was created using different numpy array values, do not use cached data
INFO:nlmod.cache:cache was created using different dictionaries, do not use cached data
INFO:nlmod.cache:caching data -> bathymetry_model_ds.nc
INFO:nlmod.mdims.mgrid:get first active modellayer for each cell in idomain
INFO:nlmod.cache:cache was created using different numpy array values, do not use cached data
INFO:nlmod.cache:cache was created using different dictionaries, do not use cached data
INFO:nlmod.mfpackages.mfpackages:creating modflow SIM, TDIS, GWF and IMS
INFO:nlmod.cache:caching data -> surface_water.nc```
Could this be that there are NaN values in the array and a different way of comparing two arrays should be used?
Is an optimization used to estimate certain values, so thtat there is always a small random error? Could we use np.isclose here?
Can't we easily make a test of this? Run a notebook twice and make sure nothing is downloaded the second time from the internet..
Best regards,
Bas
Flopy uses botm, let's stay close to the flopy naming convention.
hopefully make testing a bit more efficient?
Er zijn een aantal functies in de mgrid module met namen als data_array_3d_to_rec_list
en data_array_2d_vertex_to_rec_list
waarvan het vrij ondoorzichtig is wat ze doen en waarom ze zo'n lange naam hebben. Kijken of de structuur van deze functies verbeterd kan worden en de naamgeving logischer kan. Evt. met uitleg notebook.
Als algmene regel kan worden afgekort:
Seaborn package is not in the setup.py install requirements, nor is it in the requirements.txt
If i pass delr=delc=50.0, I get a flopy error about delr and delc not being arrays or something along those lines... I'll add the error here later. But we might need to take a look at that.
I would like a function to create a simulation object (sim
) from a model_ds
object. The ds_model is easy to store and backup, while the gwf/sim is needed to use the plot functions of flopy. I can also imagine plenty of usecases for the reversed. Probably a lot of code is already part of the package. Thanks!
The tests in github actions that try to run a modflow model encounter the following error:
mf6: /lib/x86_64-linux-gnu/libm.so.6: version 'GLIBC_2.29' not found (required by mf6)
This has probably something to do with compiling the mf6 executable in the linux test environment. I try to use the same steps as flopy does in there automated testing described here:
https://github.com/modflowpy/flopy/blob/develop/.github/workflows/ci.yml
Hi All,
It looks like the install requirements are not properly configured. Only the packages listed in setup.py > install_requires are installed when installing nlmod via pip install -e .
. The list of packages shown in Requirements.txt
is much longer. Geopandas is one of those packages that is not installed with pip install -e .
.
For a possible solution:
https://stackoverflow.com/questions/14399534/reference-requirements-txt-for-the-install-requires-kwarg-in-setuptools-setup-py
Upon configuring flopy packages (several functions within nlmod/mfpackages) check if boundary condition is placed in inactive idomain to prevent errors such as placing drains/lakes in inactive cells.
It would be nice if the DataArrays that are part of model_ds would have the attribute units
. The units are then automagically used in the labels when plotting the DataArray.
Add to the example notebook before import nlmode:
import sys
module_path = os.path.abspath(os.path.join('..'))
if module_path not in sys.path:
sys.path.append(module_path)
to make sure "import nlmod" works in the current folder structure.
In the notebook 02_surface_water the river stress period data is built from a geodataframe with:
riv_spd = nlmod.mfpackages.surface_water.build_spd(celldata, "RIV", model_ds)
The logs show that there are a lot of nan values in the rbot data. It does not result in an error but I was wondering if this is expected behavior. Also if this is the expected behavior should we remove the warning? Maybe we can only make a warning if some cells were removed because of missing/erroneous data.
The plotting functions of nlmod require vertices. It would be nice if those were added to model_ds during grid discretisations (via get_vertices
)
The mfpackages module is only meant to create modflow packages from the model dataset. Therefore the chd_at_model_edge_from_model_ds
should be split into a part that creates the 'chd' data variable for model dataset and the part that create the chd mf6 package.
Op dit moment worden de download en vergrid functies vaak in één keer aangeroepen. Bijv. in de functies nlmod.read.ahn.get_ahn()
en nlmod.read.knmi.get_recharge()
. Het zou handig zijn om deze functies van elkaar te kunnen scheiden. Dus dat je in het project eenmaal de data die je nodig hebt download en daarna naar deze data verwijst en deze vergrid bij het maken van een specifiek model.
Je kan dan bijv. het AHN van heel je projectgebied op fijne schaal downloaden (eenmalig en langzaam) en deze voor verschillende modellen vergridden (snel).
The function get_ahn_within_extent
writes a tif file to a temporary directory. It would be better to write to the cache directory. This is complex because the cache directory is only accessed via the cache decorator.
Some (suboptimal) ideas to solve this:
get_ahn
function with a cachedir_ahnI added documentation, available via rtd. There are quite some improvements to be made:
I thought this was the case before. Please have a look at the dtscalibration package if you experience any issues with the implementation.
topbot = np.concatenate([model_ds.top.data[np.newaxis], model_ds.bot.data], axis=0)
thickness = -1 * np.diff(topbot, axis=0)
idomain = xr.where(thickness < 1e-5, -1, model_ds["idomain"])
The function to combine regis and geotop add_geotop_to_regis_hlc
is very rigid. Is assumes there are no layers above the Holoceen in Regis. This usually works but I think a general function in which you can choose more layers that you want to replace with geotop is better.
Many functions have model_ds
(a model dataset) as an argument. In these functions data arrays and attributes from the model dataset are used. When you call this function it is very unclear which data arrays or attributes are necessary to run the function. We should change these functions to include all the data variables and/or attributes as function arguments. Then each function can be called using model_ds or by providing all the other function arguments.
Current function ghb_from_model_ds
in nlmod.mfpackages.mfpackages
def ghb_from_model_ds(model_ds, gwf, da_name):
"""get general head boundary from model dataset.
Parameters
----------
model_ds : xarray.Dataset
dataset with model data.
gwf : flopy ModflowGwf
groundwaterflow object.
da_name : str
name of the ghb files in the model dataset.
Proposal for new function
def ghb_from_model_ds(gwf,
model_ds=None,
ghb_cond=None,
ghb_stage=None,
gridtype=None,
only_first_active_layer=True,
first_active_layer=None,
only_active_cells=False,
idomain=None):
"""get general head boundary from model dataset.
Parameters
----------
gwf : flopy ModflowGwf
groundwaterflow object.
model_ds : xarray.Dataset or None, optional
dataset with model data. If model_ds is None the function arguments
ghb_cond, ghb_peil and gridtype should be defined. Default is None
ghb_cond : xarray.DataArray or None, optional
dataarray with the conductivity of the ghb. All non-zero cells are used
as ghb cells. If None the function argument model_ds should be defined.
Default is None
ghb_stage : xarray.DataArray or None, optional
dataarray with the stage of the ghb. If None the function argument
model_ds should be defined. Default is None
gridtype : str or None, optional
can be 'structured' or 'vertex'. If None the function argument model_ds
should be defined. Default is None
only_first_active_layer : bool, optional
if True only the ghb is only added to the first active layer in the
model. Default is True
first_active_layer : xarray.DataArray or None, optional
dataarray with for each cell the first active layer. Only necesarry
when use_first_active_layer is True. Default is None
only_active_cells : bool, optional
if True ghb cells are only to cells with an idomain of 1. Default is
False.
idomain : xarray.DataArray or None, optional
dataarray with the idomain. Only necesarry
when only_active_cells is True. Default is None
I think it is one of the more recent functions in flopy. So that the code of the get_heads_array
and get_heads_dataarray
can be written a lot leaner:
ml = sim.get_model(sim.name)
# collect heads
heads = ml.output.head().get_alldata()[:, :, 0, :]
heads[heads == ml.hdry] = np.nan
heads[heads == ml.hnoflo] = np.nan
model_ds['heads'] = xr.DataArray(
data=heads,
dims=('time', 'layer', 'cid'),
attrs={'units': 'mNAP'})
where ml
is the grondwater flow model, which is in other examples denoted with gwf
. For simple cases this also works for reading the budgetfile, but for example not for the Bergen case..
The units are also stored in the flopy objects:
length_unit = ml.disv.length_units.get_data()
time_unit = sim.tdis.time_units.get_data()
It would be nice to have all all dataarrays with units. Especially with the budgetfile this is important.
When I run the get_exes.py
script on windows with pymake = 1.2.1 I get the following error:
Updating MODFLOW 6 executables from the nightly-build repo
Attempting to download the file:
MODFLOW-USGS/modflow6-nightly-build/win64.zip
open request attempt 1 of 10
open request attempt 2 of 10
open request attempt 3 of 10
open request attempt 4 of 10
open request attempt 5 of 10
open request attempt 6 of 10
open request attempt 7 of 10
open request attempt 8 of 10
open request attempt 9 of 10
open request attempt 10 of 10
Cannot open request from:
MODFLOW-USGS/modflow6-nightly-build/win64.zip
Traceback (most recent call last):
File "C:\Anaconda3\envs\dev\lib\site-packages\pymake\utils\download.py", line 208, in _request_get
req = requests.get(
File "C:\Anaconda3\envs\dev\lib\site-packages\requests\api.py", line 76, in get
return request('get', url, params=params, **kwargs)
File "C:\Anaconda3\envs\dev\lib\site-packages\requests\api.py", line 61, in request
return session.request(method=method, url=url, **kwargs)
File "C:\Anaconda3\envs\dev\lib\site-packages\requests\sessions.py", line 516, in request
prep = self.prepare_request(req)
File "C:\Anaconda3\envs\dev\lib\site-packages\requests\sessions.py", line 449, in prepare_request
p.prepare(
File "C:\Anaconda3\envs\dev\lib\site-packages\requests\models.py", line 314, in prepare
self.prepare_url(url, params)
File "C:\Anaconda3\envs\dev\lib\site-packages\requests\models.py", line 388, in prepare_url
raise MissingSchema(error)
MissingSchema: Invalid URL 'MODFLOW-USGS/modflow6-nightly-build/win64.zip': No schema supplied. Perhaps you meant http://MODFLOW-USGS/modflow6-nightly-build/win64.zip?
During handling of the above exception, another exception occurred:
Traceback (most recent call last):
File "C:\Users\onno_\02_git_repos\nlmod\tests\get_exes.py", line 160, in <module>
test_download_nightly_build()
File "C:\Users\onno_\02_git_repos\nlmod\tests\get_exes.py", line 131, in test_download_nightly_build
pymake.download_and_unzip(url, exe_pth, verbose=True)
File "C:\Anaconda3\envs\dev\lib\site-packages\pymake\utils\download.py", line 320, in download_and_unzip
req = _request_get(
File "C:\Anaconda3\envs\dev\lib\site-packages\pymake\utils\download.py", line 218, in _request_get
req.raise_for_status()
UnboundLocalError: local variable 'req' referenced before assignment
However the test_download_and_unzip
function did run and the executables were installed on my machine. Therefore I was wondering if we need the nightly build function or we can just skip it?
Also the list_exes()
function in the get_exes.py
is only made for linux and fails in windows. I think we can remove this function as well?
In Jupyter Notebooks all logs are highlighted in red. This make the logs look like errors or warnings while often they are just information about what is happening. It would be nice if we can highlight the errors based on the log level. It might be possible to do this using the coloredlogs
package -> https://coloredlogs.readthedocs.io/en/latest/
There is already documentation available but not a readthedocs page yet.
Now you first create a model_ds_time object, and then use that to initialize a MFSimulation.
More logical to separate those steps somehow?
model_ds = nlmod.get_empty_model_ds()
followed by:
model_ds = nlmod.set_model_ds_time(...)
I'll commit something along these lines today.
In de functies regis.get_regis_dataset
en geotop.get_geotop_dataset
kan je een willekeurig extent opgeven. Echter bij sommige extents gaat het helemaal mis en krijg je in sommige cellen wel een bot maar geen kh of iets anders.
Dit is op te lossen door de extent eerst te schalen naar regis
. Dit kan met de functie regis.fit_extent_to_regis
. Maar dit moet je wel handmatig doen. Als je dat niet doet krijg je foutmeldingen. Voor nu heb ik een check ingebouwd in de eerder genoemde functies maar dit is een beetje houtje/touwtje.
Het mooiste zou zijn als je uiteindelijk ieder willekeurig extent kan meegeven in de functies regis.get_regis_dataset
en geotop.get_geotop_dataset
.
It seems easy but the 'units' attribute of the time coordinate seams to be locked and I get this error when I try to write it to a netcdf file: ValueError: failed to prevent overwriting existing key units in attrs on variable 'time'. This is probably an encoding field used by xarray to describe how a variable is serialized. To proceed, remove this key from the variable's attributes manually.
It is also mentioned here: pydata/xarray#3739. I don't have time now to dive into this. so I will revert the commit for now and create an issue for this.
Originally posted by @OnnoEbbens in #58 (comment)
For convenience:
nlmod.utils.clear_cache(model_ws)
if the keyword argument data_variables
in the function dataarray_to_shapefile
is a str an unclear error is raised. Fix function for str input
Some stuff that came up and I want to fix later:
add_knmi_to_model_dataset
to get_recharge
.find_sea_cells
from jarkus module to rws moduleget_thickness_from_topbot
and calculate_thickness
. Remove one and make more uniformda = xr.zeros_like(model_ds['top'])
. This is pretty ugly as it assumes that top is a 2d variable and all the attributes are copied to the new dataarray.If you call the util.get_cache_netcdf
function with a function (as the get_dataset_func
argument) which in turn calls the util.get_cache_netcdf
function unexpected errors may occur. The cachedir
argument that was passed to the original util.get_cache_netcdf
function is not passed onto the get_dataset_func
function. Thus when util.get_cache_netcdf
is called for a second time it is assumed that cachedir=None
. If cachedir is None
an attempt is made to create a temporary cache directory. Creating or accessing a temporary cache directory is not always possible and can result in an error. Also it is unexpected that a temporary directory is used for caching while a cache directory was specified in the original call of util.get_cache_netcdf
.
some thoughts:
util.get_cache_netcdf
function with a function (as the get_dataset_func
argument) that calls the util.get_cache_netcdf
. Why cache the same thing twice? However it is also hard to check that this occurs.cachedir=None
when calling the util.get_cache_netcdf
function. A solution will be to raise an error when cachedir=None
.Some stuff that I still want to do:
Should we make extrapolate_regis part of get_combined_layer_models? And add a keyword argument to optionally turn it off?
It sounds like a sane default to have it part of get_combined_layer_model and makes it so that users have to remember one line less.
Thanks!
The function nlmod.mlayers.update_model_ds_from_ml_layer_ds
has a lot of steps that are quite unclear. Especially the add_northsea
part should be a separate function.
the package DataSetCrossSection from art_tools is not available through pip or conda. It looks like a tool from Deltares. Please provide reference on how to access and install the package.
Using the same definitions as flopy makes sense.
It would be nice to give additional colorbar kwargs to nlmod.visualise.plots.plot_vertex_array()
. For example to shrink the colorbar.
The default directory to save mfbinaries when calling nlmod.util.download_mfbinaries()
is the parent directory of the nlmod package. This makes sense when you have the nlmod github repo because the parent directory of the nlmod package is still part of the repo. When you install nlmod using pip, the parent directory of the nlmod package is the site-packages directory. It is not a good idea to store the mfbinaries there. Therefore we should add the mfbinaries to the nlmod package directory instead.
I would be interested in an additional function in the util.py file or maybe it can be an extension of write_and_run_model
. From the OC-package we know whether head and budget are stored. The function would read and store the outputs of the model in model_ds
.
Additionally, I would like to see a option in the write_and_run_model
function to run the model in a temporary directory. Thus within
a with statement.
with tempfile.TemporaryDirectory() as tmpdirname:
# alter model_ws to tmpdirname
# write packages
# run model
# read and store outputs in model_ds using info from the OC-package
# revert model_ws
https://ugrid-conventions.github.io/ugrid-conventions/
inspired by imod-python
In de functie nlmod.read.rws.get_surface_water
wordt een shapefile ingelezen met grote oppervlaktewaterlichamen. Bij het aanroepen van de functie is niet duidelijk waar deze data vandaan komt. Zorg dat bij het aanroepen van de functie expliciet wordt vermeld welke data wordt ingelezen. Geef ook de mogelijkheid om bijv. zelf een shapefile op te geven.
Dit geldt bijvoorbeeld ook voor de conversietabel voor geotop die gebruikt wordt.
A declarative, efficient, and flexible JavaScript library for building user interfaces.
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google ❤️ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.