{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# River Profile [Py]\n", "\n", "A simple example simulating the evolution of a river longitudinal profile under the action of block-uplift vs. bedrock channel erosion, starting from a very gentle channel slope. See also {doc}`river_profile_cpp`.\n", "\n", "The local rate of elevation change, {math}`\\partial h/\\partial t`, is determined by the balance between uplift (uniform in space and time) {math}`U` and bedrock channel erosion modeled using the Stream-Power Law (SPL).\n", "\n", "```{math}\n", ":label: eq_river\n", "\\frac{\\partial h}{\\partial t} = U - K_f A^m \\left(\\frac{\\partial h}{\\partial x}\\right)^n\n", "```\n", "\n", "where {math}`A` is the drainage area and {math}`K_f`, {math}`m`, {math}`n` are SPL parameters (uniform in space and time)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import fastscapelib as fs\n", "import numpy as np\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "tags": ["remove-cell"] }, "outputs": [], "source": [ "# Theme that looks reasonably fine on both dark/light modes\n", "import matplotlib\n", "matplotlib.style.use('Solarize_Light2')\n", "matplotlib.rcParams['axes.grid'] = False" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setup the Grid, Flow Graph and Eroders\n", "\n", "Create a {py:class}`~fastscapelib.ProfileGrid` of 101 nodes each uniformly spaced by 300 meters.\n", "\n", "Set fixed-value boundary conditions on the left side and free (core) conditions on the right side of the profile." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["grid = fs.ProfileGrid(101, 300.0, [fs.NodeStatus.FIXED_VALUE, fs.NodeStatus.CORE])"] }, { "cell_type": "markdown", "metadata": {}, "source": ["Set x-coordinate values along the profile grid, which correspond to the distance from the ridge top. Start at `x0 > 0` so that the channel head effectively starts at some distance away from the ridge top."] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "x0 = 300.0\n", "length = (grid.size - 1) * grid.spacing\n", "x = np.linspace(x0 + length, x0, num=grid.size)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Create a {py:class}`~fastscapelib.FlowGraph` object with single direction flow routing (trivial case for a 1-dimensional profile).\n", "\n", "By default, base level nodes are set from fixed value boundary conditions (left node in this example)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["flow_graph = fs.FlowGraph(grid, [fs.SingleFlowRouter()])"] }, { "cell_type": "markdown", "metadata": {}, "source": ["Compute drainage area along the profile using [Hack's Law](https://en.wikipedia.org/wiki/Hack%27s_law)."] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "hack_coef = 6.69\n", "hack_exp = 1.67\n", "\n", "drainage_area = hack_coef * x**hack_exp" ] }, { "cell_type": "markdown", "metadata": {}, "source": ["Setup the SPL eroder class with a given set of parameter values."] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "spl_eroder = fs.SPLEroder(\n", " flow_graph,\n", " k_coef=1e-4,\n", " area_exp=0.5,\n", " slope_exp=1,\n", " tolerance=1e-5,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": ["## Setup Initial Conditions and External Forcing"] }, { "cell_type": "markdown", "metadata": {}, "source": ["Create an initial river profile with a very gentle, uniform slope."] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["init_elevation = (length + x0 - x) * 1e-4\n", "elevation = init_elevation"] }, { "cell_type": "markdown", "metadata": {}, "source": ["Set upflit rate as uniform (fixed value) along the profile except at the base level node (zero)."] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["uplift_rate = np.full_like(x, 1e-3)\n", "uplift_rate[0] = 0."] }, { "cell_type": "markdown", "metadata": {}, "source": ["## Run the Model"] }, { "cell_type": "markdown", "metadata": {}, "source": ["Run the model for a few thousands of time steps. At the middle of the simulation, change the value of the SPL coefficient to simulate a change in erodibility (e.g., climate change). This should create a knicpoint in the resulting river profile."] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "dt = 1e2\n", "nsteps = 4000\n", "\n", "# flow routing is required for solving SPL, although in this example\n", "# it can be computed out of the time step loop since the\n", "# flow path (single river) will not change during the simulation\n", "# (not the case, e.g., if there were closed depressions!)\n", "flow_graph.update_routes(elevation)\n", "\n", "\n", "for step in range(nsteps):\n", " # uplift (no uplift at fixed elevation boundaries)\n", " uplifted_elevation = elevation + dt * uplift_rate\n", " \n", " # abrubpt change in erodibility\n", " if step == 3000:\n", " spl_eroder.k_coef /= 4\n", "\n", " # apply channel erosion\n", " spl_erosion = spl_eroder.erode(uplifted_elevation, drainage_area, dt)\n", " \n", " # update topography\n", " elevation = uplifted_elevation - spl_erosion\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": ["## Plot the River Profile"] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fig, ax = plt.subplots(figsize=(12, 6))\n", "ax.plot(x, elevation)\n", "plt.setp(ax, xlabel=\"x\", ylabel=\"elevation\", xlim=(x[0] + 300.0, 0));" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.9.10" } }, "nbformat": 4, "nbformat_minor": 4 }