diff --git a/.github/workflows/build_docs.yaml b/.github/workflows/build_docs.yaml index d6607a5d..8766521b 100644 --- a/.github/workflows/build_docs.yaml +++ b/.github/workflows/build_docs.yaml @@ -39,6 +39,7 @@ jobs: with: path: ~/.local/bin github_token: ${{ secrets.GITHUB_TOKEN }} + repo: modflow6 - name: Modflow version shell: bash -l {0} run: mf6 --version diff --git a/.github/workflows/test.yaml b/.github/workflows/test.yaml index 288f74bd..0c30c3c3 100644 --- a/.github/workflows/test.yaml +++ b/.github/workflows/test.yaml @@ -11,13 +11,15 @@ on: jobs: build: - name: ${{ matrix.python-version }}, ${{ matrix.os }} + name: ${{ matrix.env }}, ${{ matrix.os }} runs-on: ${{ matrix.os }} strategy: fail-fast: false matrix: os: [ ubuntu-latest, macos-latest, windows-latest ] - python-version: ['3.12'] + env: + - installation/environment.yml + - ci/pandas2_env.yml steps: - name: Checkout source @@ -31,7 +33,7 @@ jobs: - name: Setup Micromamba uses: mamba-org/setup-micromamba@main with: - environment-file: installation/environment.yml + environment-file: ${{ matrix.env }} cache-environment: false cache-downloads: false # persist on the same day. @@ -51,6 +53,8 @@ jobs: with: path: ~/.local/bin github_token: ${{ secrets.GITHUB_TOKEN }} + repo: modflow6 + cache: false - name: Modflow version shell: bash -l {0} run: mf6 --version @@ -75,6 +79,6 @@ jobs: - name: Run tests shell: bash -l {0} run: | - pytest tests/test_notebooks.py -v + pytest tests/test_notebooks.py env: GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} diff --git a/AGENDA.md b/AGENDA.md index c99b23f1..8500758c 100644 --- a/AGENDA.md +++ b/AGENDA.md @@ -57,6 +57,6 @@ Monday - 1500-1515: break - 1515-1545: 08 Modflow-setup demo (Andy) - 1545-1600: break -- 1600-1630: 09 Run and process MF6T (Kalle) +- 1600-1630: 09 Run and process MF6 GWT (Kalle) - 1630-1700: 10 Run and process MODPATH/MF6-PRT models (Andy) - 1700: Wrap up diff --git a/ci/pandas2_env.yml b/ci/pandas2_env.yml new file mode 100644 index 00000000..f78dba70 --- /dev/null +++ b/ci/pandas2_env.yml @@ -0,0 +1,52 @@ +name: pyclass +channels: + - conda-forge + - nodefaults +dependencies: + # required + - python=3.14 + - pip + + - numpy + - matplotlib + + - jupyter + - jupytext + - jupyterlab + - git + - python-dateutil + - affine + - scipy + - openpyxl + - xlrd + - pandas<3 + - netcdf4 + - pyshp + - rasterio + - rasterstats + - fiona + - descartes + - pyproj + - shapely + - geos + - geojson + - geopandas + - xarray + - rioxarray + - uxarray + - pyyaml + - rtree + - pyvista + - vtk + - imageio + - requests + - pytest + - statsmodels + - dataretrieval + - flopy=3.9 + - gis-utils + - sfrmaker + - contextily + - cartopy + - modflow-export + - modflow-setup diff --git a/docs/docs-environment.yml b/docs/docs-environment.yml index 35a034d6..c653f7e3 100644 --- a/docs/docs-environment.yml +++ b/docs/docs-environment.yml @@ -3,7 +3,7 @@ channels: - conda-forge dependencies: # required - - python=3.11 + - python=3.14 - pip - numpy diff --git a/installation/environment.yml b/installation/environment.yml index fefb3afb..3b36cc43 100644 --- a/installation/environment.yml +++ b/installation/environment.yml @@ -19,7 +19,7 @@ dependencies: - scipy - openpyxl - xlrd - - pandas<3.0 + - pandas - netcdf4 - pyshp - rasterio @@ -43,7 +43,7 @@ dependencies: - pytest - statsmodels - dataretrieval - - flopy>=3.7 + - flopy=3.9 - gis-utils - sfrmaker - contextily diff --git a/notebooks/part1_flopy/03_Loading_and_visualizing_models.ipynb b/notebooks/part1_flopy/03_Loading_and_visualizing_models.ipynb index 24dbfcf4..fc4e7807 100644 --- a/notebooks/part1_flopy/03_Loading_and_visualizing_models.ipynb +++ b/notebooks/part1_flopy/03_Loading_and_visualizing_models.ipynb @@ -183,22 +183,6 @@ "ax.set_ylabel(f'{modelgrid.units.capitalize()} northing, {modelgrid.crs.name}')" ] }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Exporting the model grid to a shapefile" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "modelgrid.write_shapefile(str(output_folder / 'pleasant_grid.shp'))" - ] - }, { "cell_type": "markdown", "metadata": {}, @@ -213,41 +197,15 @@ "outputs": [], "source": [ "modelgrid = model.modelgrid\n", - "modelgrid_gdf = modelgrid.geo_dataframe" + "modelgrid_gdf = modelgrid.geo_dataframe\n", + "modelgrid_gdf" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "#### Add row, column information to model grid `GeoDataFrame`" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "modelgrid.shape" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "modelgrid_gdf.index[-1] // modelgrid.ncol" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "1 // modelgrid.ncol" + "### Exporting the model grid to a shapefile" ] }, { @@ -256,9 +214,7 @@ "metadata": {}, "outputs": [], "source": [ - "modelgrid_gdf['i'] = modelgrid_gdf.index // modelgrid.ncol\n", - "modelgrid_gdf['j'] = modelgrid_gdf.index % modelgrid.ncol\n", - "modelgrid_gdf.head()" + "modelgrid_gdf.to_file(output_folder / 'pleasant_grid.shp')" ] }, { @@ -1100,6 +1056,10 @@ "fig, ax = plt.subplots(figsize=(10, 5))\n", "in_cols = ['STO-SS_IN', 'STO-SY_IN', 'WEL_IN', 'RCHA_IN', 'CHD_IN', 'SFR_IN', 'LAK_IN']\n", "out_cols = [c.replace('_IN', '_OUT') for c in in_cols]\n", + "# as of pandas 3.0.1, \n", + "# the datatime index of the flux dataframe produced by Flopy produces an error with .plot.bar();\n", + "# switch the index to a categorical (date string) index instead\n", + "flux.index = flux.index.strftime('%Y-%m-%d')\n", "flux[in_cols].plot.bar(stacked=True, ax=ax)\n", "(-flux[out_cols]).plot.bar(stacked=True, ax=ax)\n", "ax.legend(loc='lower left', bbox_to_anchor=(1, 0))\n", @@ -1107,6 +1067,13 @@ "ax.set_ylabel('Simulated Flux, in $m^3/d$')" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, { "cell_type": "markdown", "metadata": {}, @@ -1142,7 +1109,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.11" + "version": "3.14.3" } }, "nbformat": 4, diff --git a/notebooks/part1_flopy/04_Modelgrid_and_intersection.ipynb b/notebooks/part1_flopy/04_Modelgrid_and_intersection.ipynb index 50f5bc28..34554ee3 100644 --- a/notebooks/part1_flopy/04_Modelgrid_and_intersection.ipynb +++ b/notebooks/part1_flopy/04_Modelgrid_and_intersection.ipynb @@ -1202,7 +1202,7 @@ "outputs": [], "source": [ "main_stem = []\n", - "gx = GridIntersect(modelgrid, method=\"vertex\")\n", + "gx = GridIntersect(modelgrid)\n", "for iloc, row in flw_main.iterrows():\n", " results = gx.intersect(row.geometry)\n", " main_stem.extend(results.cellids)\n", diff --git a/notebooks/part1_flopy/05_Unstructured_Grid_generation.ipynb b/notebooks/part1_flopy/05_Unstructured_Grid_generation.ipynb index 0c4d76ac..bd4decd5 100644 --- a/notebooks/part1_flopy/05_Unstructured_Grid_generation.ipynb +++ b/notebooks/part1_flopy/05_Unstructured_Grid_generation.ipynb @@ -399,7 +399,7 @@ "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(sgrid, method=\"vertex\")\n", + "ix = GridIntersect(sgrid)\n", "result = ix.intersect(basindf.geometry.values[0], contains_centroid=True)\n", "rowcol = result[\"cellids\"]" ] @@ -550,7 +550,7 @@ "outputs": [], "source": [ "rowcols = []\n", - "ix = GridIntersect(sgrid, method=\"vertex\")\n", + "ix = GridIntersect(sgrid)\n", "for geom in nhddf.geometry.values:\n", " rcs = ix.intersects(geom)[\"cellids\"]\n", " rowcols.extend(list(rcs))" @@ -673,7 +673,7 @@ "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(rfgrid, method=\"vertex\")\n", + "ix = GridIntersect(rfgrid)\n", "result = ix.intersect(basindf.geometry.values[0], contains_centroid=True)\n", "rowcol = result[\"cellids\"]\n", "\n", @@ -849,7 +849,7 @@ "# live code this/give it as an activity\n", "idomainp = np.ones(sgrid.idomain.shape, dtype=int)\n", "\n", - "gix = GridIntersect(sgrid, method=\"vertex\")\n", + "gix = GridIntersect(sgrid)\n", "result = gix.intersect(lgr_gdf.geometry.values[0])\n", "cid = result[\"cellids\"]\n", "\n", @@ -959,7 +959,7 @@ "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(childgrid, method=\"vertex\")\n", + "ix = GridIntersect(childgrid)\n", "cstrmcells = []\n", "for geom in nhddf.geometry.values:\n", " rcs = ix.intersects(geom)[\"cellids\"]\n", @@ -1120,7 +1120,7 @@ "source": [ "# create the idomain array\n", "idomain = np.zeros(quadgrid.shape, dtype=int)\n", - "ix = GridIntersect(quadgrid, method=\"vertex\")\n", + "ix = GridIntersect(quadgrid)\n", "nodes = ix.intersect(basindf.geometry.values[0], contains_centroid=True)[\"cellids\"]\n", "\n", "idomain[:, list(nodes)] = 1\n", @@ -1136,7 +1136,7 @@ "source": [ "# create an array of stream cells\n", "qstr = []\n", - "ix = GridIntersect(quadgrid, method=\"vertex\")\n", + "ix = GridIntersect(quadgrid)\n", "for geom in nhddf.geometry.values:\n", " nodes = ix.intersects(geom)[\"cellids\"]\n", " qstr.extend(list(nodes))" @@ -1383,7 +1383,7 @@ "outputs": [], "source": [ "tristr = []\n", - "ix = GridIntersect(trigrid, method=\"vertex\")\n", + "ix = GridIntersect(trigrid)\n", "for geom in nhddf.geometry.values:\n", " nodes = ix.intersects(geom)[\"cellids\"]\n", " tristr.extend(list(nodes))\n", @@ -1519,7 +1519,7 @@ "outputs": [], "source": [ "vorst = []\n", - "ix = GridIntersect(vorgrid, method=\"vertex\")\n", + "ix = GridIntersect(vorgrid)\n", "for geom in nhddf.geometry.values:\n", " nodes = ix.intersects(geom)[\"cellids\"]\n", " vorst.extend(list(nodes))\n", diff --git a/notebooks/part1_flopy/06-Project-structured.ipynb b/notebooks/part1_flopy/06-Project-structured.ipynb index 0d562674..7c22a148 100644 --- a/notebooks/part1_flopy/06-Project-structured.ipynb +++ b/notebooks/part1_flopy/06-Project-structured.ipynb @@ -143,7 +143,7 @@ "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(base_grid, method=\"vertex\", rtree=True)" + "ix = GridIntersect(base_grid, rtree=True)" ] }, { diff --git a/notebooks/part1_flopy/07-stream_capture_voronoi.ipynb b/notebooks/part1_flopy/07-stream_capture_voronoi.ipynb index bddf6030..cbea80d7 100644 --- a/notebooks/part1_flopy/07-stream_capture_voronoi.ipynb +++ b/notebooks/part1_flopy/07-stream_capture_voronoi.ipynb @@ -661,8 +661,8 @@ "for gage_num in range(2):\n", " for lay_num in range(3):\n", " cax = ax[ax_num]\n", - " data = drf.loc[lay_num][f'Gage{gage_num+1}'].values\n", - " data[data<0] =0\n", + " data = drf.loc[lay_num][f'Gage{gage_num+1}'].values.copy()\n", + " data[data < 0] = 0\n", " mm = fp.plot.PlotMapView(model= base_m, ax=cax)\n", " mm.plot_bc('SFR', plotAll=True, color='blue')\n", " mm.plot_bc('WEL', plotAll=True)\n", @@ -696,9 +696,21 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3 (ipykernel)", + "display_name": "pyclass", "language": "python", "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.14.3" } }, "nbformat": 4, diff --git a/notebooks/part1_flopy/10a_prt_particle_tracking-demo.ipynb b/notebooks/part1_flopy/10a_prt_particle_tracking-demo.ipynb index 8ec23325..b646ff37 100644 --- a/notebooks/part1_flopy/10a_prt_particle_tracking-demo.ipynb +++ b/notebooks/part1_flopy/10a_prt_particle_tracking-demo.ipynb @@ -133,7 +133,8 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "### Change the workspace" + "### Change the workspace\n", + "(so that we don't modify the original flow model)" ] }, { @@ -180,6 +181,172 @@ "sim.write_simulation()" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Add `IFLOWFACE=6` information to the SFR Package\n", + "This may not be necessary here because of the quadtree grid (that is refined in the stream cells), but it is often best practice.\n", + "\n", + "Particle tracking with numerical models is prone to the problem of weak sinks, where a sink (flux of water leaving a cell) is not strong enough to capture all of the flow through that cell. In other words, one or more cell faces still has outward flow, allowing water to pass through the cell. By default, PRT and MODPATH assume that sinks are distributed uniformly throughout a cell (`IFLOWFACE=0` or `IFACE=0`; see below). This setting is appropriate for wells, but for surface water features (especially small streams), this often prevents particles from discharging, leading to unrealistic pathlines and travel times. Some possible solutions:\n", + "* reduce the grid size until all cells of interest are strong sinks\n", + "* set `STOP_AT_WEAK_SINK=True` in the PRT Particle Release (PRP) input (or the equivalent setting in MODPATH), causing all particles to terminate at weak sinks. This option is often not appropriate for small streams, which in many cases should act as weak sinks (only capturing the shallowest particles).\\\n", + "* set `IFLOWFACE=-1` or `IFACE=6` for cells with surface water boundary conditions, to instruct PRT or MODPATH to apply the sink discharge to the top face. In cells with small (weak) streams, assuming the stream bottom corresponds to the cell top, shallow particles that are near the top of the cell will discharge, while deeper particles will be allowed to pass through realistically. Abrams and others (2013) discuss this concept in more detail.\n", + "\n", + "**From the MODFLOW 6 I/O guide:**\n", + "By default, flows associated with stress packages are assumed to be distributed uniformly over the volume of a cell.\n", + "Distributed external inflows and outflows are reflected in the cell-cell flows calculated by the GWF Model, but they are\n", + "not directly involved in determining the normal face velocities used to track a particle through the cell. The user can\n", + "optionally assign a flow associated with a stress package to any face of the cell. Assignment of external flows is done by\n", + "setting the value of an auxiliary package variable called IFLOWFACE to a value that corresponds to one of the cell faces.\n", + "To accommodate non-rectangular cells, the IFLOWFACE numbering scheme in the PRT Model is based on clockwise\n", + "ordering of the lateral cell faces. For a DIS-grid cell, IFLOWFACE = 1 corresponds to the “western” face, i.e., the face parallel to the y axis that has the lesser x coordinate. For a DISV-grid cell, IFLOWFACE = 1 corresponds to the face\n", + "that extends in the clockwise direction from the first vertex listed for that cell in the CELL2D input block. IFLOWFACE\n", + "numbering of lateral cell faces proceeds in clockwise direction. IFLOWFACE = -2 corresponds to the cell bottom, and\n", + "IFLOWFACE = -1 corresponds to the cell top.\n", + "\n", + "**In MODPATH 7** (and earlier)\n", + "`IFACE` is equivalent to `IFLOWFACE` in MODFLOW 6 PRT. Default `IFACE` values can be specified by package, or `IFACE` can be specified by cell (see the docs). A key difference is that `IFACE=6` corresponds to the top of the cell, and `IFACE=5` corresponds to the bottom." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Two ways to add `IFLOWFACE`\n", + "\n", + "**1) Using Flopy**\n", + "\n", + "Add an auxiliary column named `'iflowface'` to the Options block" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from copy import deepcopy\n", + "\n", + "gwf2 = deepcopy(gwf)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "gwf2.sfr.auxiliary = ['iflowface']\n", + "gwf2.sfr.auxiliary" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Add `iflowface` to the package data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "packagedata = pd.DataFrame(gwf2.sfr.packagedata.array)\n", + "# insert the auxiliary column before the existing boundname column\n", + "# (so that the boundname column is still last)\n", + "packagedata.insert(len(packagedata.columns)-1, 'iflowface', -1)\n", + "# make sure the iflowface column is an integer!\n", + "packagedata['iflowface'] = packagedata['iflowface'].astype(int)\n", + "# cast packagedata back to a recarray (required by Flopy)\n", + "# Note: index=False is needed here, \n", + "# otherwise the extra index column inserted by pandas edges out the boundnames column\n", + "# resulting in no boundnames being passed to Flopy\n", + "gwf2.sfr.packagedata.set_data(packagedata.to_records(index=False))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "gwf2.sfr.packagedata.array[:5]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**Alternatively, we can write the model files, and then edit the SFR Package file directly.**" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "write_iflowface = True\n", + "set_iflowface = -1\n", + "#original_sfr_filename = load_ws / gwf.sfr.filename\n", + "sfr_filename = gwf_ws / gwf.sfr.filename\n", + "original_sfr_filename = sfr_filename\n", + "new_sfr_lines = list()\n", + "package_already_has_iflowface = False\n", + "aux_entry_added = False\n", + "boundnames = False\n", + "with open(original_sfr_filename) as src:\n", + " for line in src:\n", + " normalized_line = ' '.join(line.lower().split()).lower()\n", + " if 'begin options' in normalized_line:\n", + " new_sfr_lines.append(line)\n", + " for line in src:\n", + " normalized_line = ' '.join(line.lower().split()).lower()\n", + " if 'boundnames' in normalized_line:\n", + " boundnames = True\n", + " if 'auxiliary' in normalized_line:\n", + " if 'iflowface' not in normalized_line:\n", + " items = line.strip().split()\n", + " items.append('iflowface')\n", + " new_sfr_lines.append(f\" {' '.join(items)}\\n\")\n", + " else:\n", + " new_sfr_lines.append(line)\n", + " package_already_has_iflowface = True\n", + " aux_entry_added = True\n", + " elif not aux_entry_added and 'end options' in normalized_line:\n", + " new_sfr_lines.append(' auxiliary iflowface\\n')\n", + " new_sfr_lines.append(line)\n", + " aux_entry_added = True\n", + " break\n", + " elif 'end options' in normalized_line:\n", + " break\n", + " else:\n", + " new_sfr_lines.append(line)\n", + " elif 'begin packagedata' in normalized_line:\n", + " new_sfr_lines.append(line)\n", + " for line in src:\n", + " if 'end packagedata' in line.lower():\n", + " new_sfr_lines.append(line)\n", + " break\n", + " items = line.strip().split()\n", + " # if there's a boundname column, it needs to be last\n", + " if boundnames:\n", + " items.insert(-1, str(set_iflowface))\n", + " else:\n", + " items.append(set_iflowface)\n", + " new_sfr_lines.append(f\" {' '.join(items)}\\n\")\n", + " else:\n", + " new_sfr_lines.append(line)\n", + "\n", + "if write_iflowface and not package_already_has_iflowface:\n", + " with open(sfr_filename, 'w') as dest:\n", + " for line in new_sfr_lines:\n", + " dest.write(line)" + ] + }, { "cell_type": "markdown", "metadata": {}, @@ -431,9 +598,9 @@ "outputs": [], "source": [ "sd = flopy.modpath.CellDataType()\n", - "mp7_particle_data = flopy.modpath.NodeParticleData(subdivisiondata=[sd],\n", - " nodes=list(particle_start_nodes))\n", - "prt_particle_data = list(mp7_particle_data.to_prp(prt.modelgrid))\n", + "mp7_particle_data = flopy.modpath.NodeParticleData(subdivisiondata=sd,\n", + " nodes=[cell2d for cell2d in particle_start_nodes])\n", + "prt_particle_data = list(mp7_particle_data.to_prp(prt.modelgrid, localz=True))\n", "prt_particle_data" ] }, @@ -465,6 +632,7 @@ "flopy.mf6.ModflowPrtprp(\n", " prt,\n", " nreleasepts=len(prt_start_data),\n", + " #packagedata=prt_particle_data,\n", " packagedata=prt_start_data.to_records(index=False).tolist(),\n", " local_z=True,\n", " boundnames=True,\n", @@ -483,12 +651,13 @@ "\n", "**Notes:**\n", "PRT outputs a record each time a:\n", - " 0: particle was released\n", - " 1: particle exited a cell\n", - " 2: time step ended\n", - " 3: particle terminated\n", - " 4: particle entered a weak sink cell\n", - " 5: user-specified tracking time\n", + "\n", + " 0) particle was released\n", + " 1) particle exited a cell\n", + " 2) time step ended\n", + " 3) particle terminated\n", + " 4) particle entered a weak sink cell\n", + " 5) user-specified tracking time\n", " \n", "The code below shows how to input user-specified tracking times to Flopy. Depending on the problem, this may not be necessary in a real-world context, as there may already be many records from particles exiting cells. Therefore, `tracktimes=None` (no user-specified times) is ultimately input to Flopy." ] @@ -541,7 +710,7 @@ " prt_sim,\n", " filename=f\"{prt_name}.ems\",\n", ")\n", - "sim.register_solution_package(ems, [prt.name])" + "prt_sim.register_solution_package(ems, [prt.name])" ] }, { @@ -728,9 +897,9 @@ "metadata": {}, "outputs": [], "source": [ - "df['t_years'] = df['t'] * 325.25\n", - "df['k'] = df['icell'].values[-1] // ncells_per_layer\n", - "df['cell2d'] = df['icell'].values[-1] % ncells_per_layer\n" + "df['t_years'] = df['t'] / 325.25\n", + "df['k'] = df['icell'] // ncells_per_layer\n", + "df['cell2d'] = df['icell'] % ncells_per_layer\n" ] }, { @@ -811,7 +980,9 @@ "for package, izone in izones.items():\n", " layer_name = f\"Pathlines going to {package.upper()} cells\"\n", " izone_df = lines.loc[lines['izone'] == izone]\n", - " izone_df.to_file(output_geopackage, index=False, layer=layer_name)" + " izone_df.to_file(output_geopackage, index=False, layer=layer_name)\n", + "izone_df.loc[izone_df['izone'] == 0].to_file(\n", + " output_geopackage, index=False, layer='All other pathlines')" ] }, { @@ -838,9 +1009,34 @@ " izone_df = lines.loc[lines['izone'] == izone]\n", " izone_df.plot(ax=ax, ec=line_colors[package],\n", " label=f'Pathlines going to {package.upper()} cells')\n", + "lines.loc[lines['izone'] == 0].plot(\n", + " ax=ax, ec='g',\n", + " label=f'All other pathlines')\n", "pmv.plot_ibound()\n", "pmv.ax.legend()" ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "lines" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### References\n", + "Abrams, D., Haitjema, H. and Kauffman, L. (2013), On Modeling Weak Sinks in MODPATH. Groundwater, 51: 597-602. https://doi.org/10.1111/j.1745-6584.2012.00995.x" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [] } ], "metadata": { @@ -848,6 +1044,18 @@ "display_name": "pyclass", "language": "python", "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.14.3" } }, "nbformat": 4, diff --git a/notebooks/part1_flopy/solutions/03_Loading_and_visualizing_models-solutions.ipynb b/notebooks/part1_flopy/solutions/03_Loading_and_visualizing_models-solutions.ipynb index 61df9cbe..05ff75fb 100644 --- a/notebooks/part1_flopy/solutions/03_Loading_and_visualizing_models-solutions.ipynb +++ b/notebooks/part1_flopy/solutions/03_Loading_and_visualizing_models-solutions.ipynb @@ -182,22 +182,6 @@ "ax.set_ylabel(f'{modelgrid.units.capitalize()} northing, {modelgrid.crs.name}')" ] }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Exporting the model grid to a shapefile" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "modelgrid.write_shapefile(str(output_folder / 'pleasant_grid.shp'))" - ] - }, { "cell_type": "markdown", "metadata": {}, @@ -219,16 +203,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "#### Add row, column information to model grid `GeoDataFrame`" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "modelgrid.shape" + "### Exporting the model grid to a shapefile" ] }, { @@ -237,9 +212,7 @@ "metadata": {}, "outputs": [], "source": [ - "modelgrid_gdf['i'] = modelgrid_gdf.index // modelgrid.ncol\n", - "modelgrid_gdf['j'] = modelgrid_gdf.index % modelgrid.ncol\n", - "modelgrid_gdf" + "modelgrid_gdf.to_file(output_folder / 'pleasant_grid.shp')" ] }, { @@ -1021,7 +994,7 @@ "metadata": {}, "outputs": [], "source": [ - "pd.read_csv(sim_ws / 'external/chd_001.dat', delim_whitespace=True)" + "pd.read_csv(sim_ws / 'external/chd_001.dat', sep=r'\\s+')" ] }, { @@ -1130,6 +1103,10 @@ "fig, ax = plt.subplots(figsize=(10, 5))\n", "in_cols = ['STO-SS_IN', 'STO-SY_IN', 'WEL_IN', 'RCHA_IN', 'CHD_IN', 'SFR_IN', 'LAK_IN']\n", "out_cols = [c.replace('_IN', '_OUT') for c in in_cols]\n", + "# as of pandas 3.0.1, \n", + "# the datatime index of the flux dataframe produced by Flopy produces an error with .plot.bar();\n", + "# switch the index to a categorical (date string) index instead\n", + "flux.index = flux.index.strftime('%Y-%m-%d')\n", "flux[in_cols].plot.bar(stacked=True, ax=ax)\n", "(-flux[out_cols]).plot.bar(stacked=True, ax=ax)\n", "ax.legend(loc='lower left', bbox_to_anchor=(1, 0))\n", @@ -1172,7 +1149,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.11" + "version": "3.14.3" } }, "nbformat": 4, diff --git a/notebooks/part1_flopy/solutions/05_Unstructured_Grid_generation.ipynb b/notebooks/part1_flopy/solutions/05_Unstructured_Grid_generation.ipynb index 839b7ed7..59160fed 100644 --- a/notebooks/part1_flopy/solutions/05_Unstructured_Grid_generation.ipynb +++ b/notebooks/part1_flopy/solutions/05_Unstructured_Grid_generation.ipynb @@ -906,7 +906,7 @@ } ], "source": [ - "ix = GridIntersect(sgrid, method=\"vertex\")\n", + "ix = GridIntersect(sgrid)\n", "result = ix.intersect(basindf.geometry.values[0], contains_centroid=True)\n", "rowcol = result[\"cellids\"]\n", "rowcol" @@ -1149,7 +1149,7 @@ "outputs": [], "source": [ "rowcols = []\n", - "ix = GridIntersect(sgrid, method=\"vertex\")\n", + "ix = GridIntersect(sgrid)\n", "for geom in nhddf.geometry.values:\n", " rcs = ix.intersects(geom)[\"cellids\"]\n", " rowcols.extend(list(rcs))" @@ -1309,7 +1309,7 @@ "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(rfgrid, method=\"vertex\")\n", + "ix = GridIntersect(rfgrid)\n", "result = ix.intersect(basindf.geometry.values[0], contains_centroid=True)\n", "rowcol = result[\"cellids\"]\n", "\n", @@ -1552,7 +1552,7 @@ "# live code this/give it as an activity\n", "idomainp = np.ones(sgrid.idomain.shape, dtype=int)\n", "\n", - "gix = GridIntersect(sgrid, method=\"vertex\")\n", + "gix = GridIntersect(sgrid)\n", "result = gix.intersect(lgr_gdf.geometry.values[0])\n", "cid = result[\"cellids\"]\n", "\n", @@ -1682,7 +1682,7 @@ "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(childgrid, method=\"vertex\")\n", + "ix = GridIntersect(childgrid)\n", "cstrmcells = []\n", "for geom in nhddf.geometry.values:\n", " rcs = ix.intersects(geom)[\"cellids\"]\n", @@ -1865,7 +1865,7 @@ "source": [ "# create the idomain array\n", "idomain = np.zeros(quadgrid.shape, dtype=int)\n", - "ix = GridIntersect(quadgrid, method=\"vertex\")\n", + "ix = GridIntersect(quadgrid)\n", "nodes = ix.intersect(basindf.geometry.values[0], contains_centroid=True)[\"cellids\"]\n", "\n", "idomain[:, list(nodes)] = 1\n", @@ -1881,7 +1881,7 @@ "source": [ "# create an array of stream cells\n", "qstr = []\n", - "ix = GridIntersect(quadgrid, method=\"vertex\")\n", + "ix = GridIntersect(quadgrid)\n", "for geom in nhddf.geometry.values:\n", " nodes = ix.intersects(geom)[\"cellids\"]\n", " qstr.extend(list(nodes))" @@ -2215,7 +2215,7 @@ "outputs": [], "source": [ "tristr = []\n", - "ix = GridIntersect(trigrid, method=\"vertex\")\n", + "ix = GridIntersect(trigrid)\n", "for geom in nhddf.geometry.values:\n", " nodes = ix.intersects(geom)[\"cellids\"]\n", " tristr.extend(list(nodes))\n", @@ -2362,7 +2362,7 @@ "outputs": [], "source": [ "vorst = []\n", - "ix = GridIntersect(vorgrid, method=\"vertex\")\n", + "ix = GridIntersect(vorgrid)\n", "for geom in nhddf.geometry.values:\n", " nodes = ix.intersects(geom)[\"cellids\"]\n", " vorst.extend(list(nodes))\n", diff --git a/notebooks/part1_flopy/solutions/06-Project-quadtree.ipynb b/notebooks/part1_flopy/solutions/06-Project-quadtree.ipynb index 5815dbd8..58e1794d 100644 --- a/notebooks/part1_flopy/solutions/06-Project-quadtree.ipynb +++ b/notebooks/part1_flopy/solutions/06-Project-quadtree.ipynb @@ -411,7 +411,7 @@ "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(base_grid, method=\"vertex\", rtree=True)" + "ix = GridIntersect(base_grid, rtree=True)" ] }, { diff --git a/notebooks/part1_flopy/solutions/06-Project-structured_completed.ipynb b/notebooks/part1_flopy/solutions/06-Project-structured_completed.ipynb index 90fd90f3..72fd5607 100644 --- a/notebooks/part1_flopy/solutions/06-Project-structured_completed.ipynb +++ b/notebooks/part1_flopy/solutions/06-Project-structured_completed.ipynb @@ -160,12 +160,12 @@ }, { "cell_type": "code", - "execution_count": 8, + "execution_count": null, "id": "9eb9a7bd-02f1-4966-b522-4bdbabeb0d40", "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(base_grid, method=\"vertex\", rtree=True)" + "ix = GridIntersect(base_grid, rtree=True)" ] }, { diff --git a/notebooks/part1_flopy/solutions/06-Project-voronoi.ipynb b/notebooks/part1_flopy/solutions/06-Project-voronoi.ipynb index 7ea74d6a..ac98b469 100644 --- a/notebooks/part1_flopy/solutions/06-Project-voronoi.ipynb +++ b/notebooks/part1_flopy/solutions/06-Project-voronoi.ipynb @@ -305,12 +305,12 @@ }, { "cell_type": "code", - "execution_count": 14, + "execution_count": null, "id": "9eb9a7bd-02f1-4966-b522-4bdbabeb0d40", "metadata": {}, "outputs": [], "source": [ - "ix = GridIntersect(base_grid, method=\"vertex\", rtree=True)" + "ix = GridIntersect(base_grid, rtree=True)" ] }, { diff --git a/notebooks/part1_flopy/solutions/07-stream_capture_voronoi.ipynb b/notebooks/part1_flopy/solutions/07-stream_capture_voronoi.ipynb index 87f20622..a3dcd5b5 100644 --- a/notebooks/part1_flopy/solutions/07-stream_capture_voronoi.ipynb +++ b/notebooks/part1_flopy/solutions/07-stream_capture_voronoi.ipynb @@ -13158,7 +13158,7 @@ }, { "cell_type": "code", - "execution_count": 47, + "execution_count": null, "id": "fd8e62ba", "metadata": {}, "outputs": [ @@ -13180,8 +13180,8 @@ "for gage_num in range(2):\n", " for lay_num in range(3):\n", " cax = ax[ax_num]\n", - " data = drf.loc[lay_num][f'Gage{gage_num+1}'].values\n", - " data[data<0] =0\n", + " data = drf.loc[lay_num][f'Gage{gage_num+1}'].values.copy()\n", + " data[data < 0] = 0\n", " mm = fp.plot.PlotMapView(model= base_m, ax=cax)\n", " mm.plot_bc('SFR', plotAll=True, color='blue')\n", " mm.plot_bc('WEL', plotAll=True)\n",