{ "cells": [ { "cell_type": "markdown", "id": "be12cf16-dfa2-4a42-8e4e-a5bfafb12389", "metadata": {}, "source": [ "### THE DANGER OF EXCESSIVELY LARGE SINGLE TIME STEPS IN DIFFUSION\n", "\n", "-> When the time step = (0.5 / diffusion rate), \n", " a 2-bin system equilibrates in a single step, \n", " and some 3-bin systems can over-shoot equilibrium!\n", "\n", "-> When the time step = (0.33333 / diffusion rate), \n", " some 3-bin systems equilibrate in a single step\n", "\n", "So, (0.33333 / diffusion rate) is a - rather lax - upper bound for\n", "sensible single time steps! \n", "\n", "IMPORTANT: The above is for delta_x = 1; in general, multiply by delta_x**2\n", "\n", "The **\"Von Neumann stability analysis\"**, which provides a slighly looser max value of time steps, is also discussed.\n", "\n", "That value of 0.33333 is saved in the Class variable \"time_step_threshold\"" ] }, { "cell_type": "markdown", "id": "f077ac33-a00f-44b3-9f60-3751acb5b8e4", "metadata": {}, "source": [ "### TAGS : \"diffusion 1D\", \"under-the-hood\"" ] }, { "cell_type": "code", "execution_count": 1, "id": "f41dc3ae-42b8-4356-9460-99c65397c3d2", "metadata": {}, "outputs": [], "source": [ "LAST_REVISED = \"May 27, 2025\"\n", "LIFE123_VERSION = \"1.0.0rc3\" # Library version this experiment is based on" ] }, { "cell_type": "code", "execution_count": 2, "id": "fcdcc8cf-6b07-4e37-9c80-8e261e0e9def", "metadata": { "tags": [] }, "outputs": [], "source": [ "#import set_path # Using MyBinder? Uncomment this before running the next cell!" ] }, { "cell_type": "code", "execution_count": 3, "id": "fe25e919", "metadata": { "tags": [] }, "outputs": [], "source": [ "#import sys, os\n", "#os.getcwd()\n", "#sys.path.append(\"C:/some_path/my_env_or_install\") # CHANGE to the folder containing your venv or libraries installation!\n", "# NOTE: If any of the imports below can't find a module, uncomment the lines above, or try: import set_path\n", "\n", "from life123 import ChemData, BioSim1D, check_version" ] }, { "cell_type": "code", "execution_count": 4, "id": "49566418-5a40-4b6b-a0c3-2c4cc021d944", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "OK\n" ] } ], "source": [ "check_version(LIFE123_VERSION)" ] }, { "cell_type": "code", "execution_count": null, "id": "229a0ef9-0efd-4a74-9899-b053bbf32746", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "bacef70f-15ac-42f7-b2b1-d74362919af7", "metadata": {}, "source": [ "## Simulate a 2-bin system, with 1 chemical" ] }, { "cell_type": "code", "execution_count": 5, "id": "61eefdac-4369-45d4-aeb4-513ef0a9c19f", "metadata": { "tags": [] }, "outputs": [], "source": [ "chem_data = ChemData(diffusion_rates=10.)" ] }, { "cell_type": "code", "execution_count": 6, "id": "992cf517-dde9-4341-8640-492aaffcf986", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 0:\n", "2 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1
0A10.0100.00.0
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1\n", "0 A 10.0 100.0 0.0" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio = BioSim1D(n_bins=2, chem_data=chem_data)\n", "bio.inject_conc_to_bin(bin_address=0, delta_conc=100., chem_index=0)\n", "bio.describe_state()" ] }, { "cell_type": "code", "execution_count": 7, "id": "c29b3617-7913-4cf2-be84-7e8def8e4476", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 0.05:\n", "2 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1
0A10.050.050.0
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1\n", "0 A 10.0 50.0 50.0" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio.time_step_threshold = 0.51 # TO BYPASS THE SAFETY LIMIT \n", " # (for the allowed maximum delta Time) that is typically in place\n", "\n", "# When the time step is (0.5 / diffusion rate),\n", "# a 2-bin system equilibrates in a single step!\n", "bio.diffuse(time_step=0.05, n_steps=1)\n", "bio.describe_state()" ] }, { "cell_type": "markdown", "id": "4d9d808c-a64e-4e75-a173-e2f46918712b", "metadata": {}, "source": [ "_Note the [50. 50.] values : **the two bins have equilibratedin a single simulation step!**_ \n", "Mighty suspicious!!" ] }, { "cell_type": "code", "execution_count": null, "id": "39928615-dc12-462e-94b1-ed82b18ccd58", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "104fce42-1ae2-456e-8160-73f48c2040b2", "metadata": {}, "source": [ "## Start over with a 3-bin system" ] }, { "cell_type": "code", "execution_count": 8, "id": "5c1bd9a7-9323-4b8a-9130-40dba22dc2ac", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 0:\n", "3 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1Bin 2
0A10.00.0100.00.0
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1 Bin 2\n", "0 A 10.0 0.0 100.0 0.0" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio = BioSim1D(n_bins=3, chem_data=chem_data)\n", "bio.inject_conc_to_bin(bin_address=1, delta_conc=100., chem_index=0)\n", "bio.describe_state()" ] }, { "cell_type": "code", "execution_count": 9, "id": "95397550-9261-494d-8c17-b410a2474b61", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 0.05:\n", "3 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1Bin 2
0A10.050.00.050.0
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1 Bin 2\n", "0 A 10.0 50.0 0.0 50.0" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio.time_step_threshold = 0.51 # TO BYPASS THE SAFETY LIMIT \n", " # (for the allowed maximum delta Time) that is typically in place\n", "\n", "# When the time step is (0.5 / diffusion rate),\n", "# a 3-bin system can overshoot equilibrium!\n", "bio.diffuse(time_step=0.05, n_steps=1)\n", "bio.describe_state()" ] }, { "cell_type": "markdown", "id": "718a6239-aa63-4574-9ef2-2b06a934844c", "metadata": {}, "source": [ "_Note the concentrations of [50. 0. 50.] : the diffusion has **over-shot equilibrium!!!**_\n", "#### Cleary, a very excessively large time step!" ] }, { "cell_type": "code", "execution_count": null, "id": "cb658630-baae-4f35-945d-2f377e980d32", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "bab05738-3e8d-4c8c-9c5a-cad9dc26ac68", "metadata": {}, "source": [ "## Start over again with a new 3-bin system" ] }, { "cell_type": "code", "execution_count": 10, "id": "3ebadd9c-3604-43d0-abe8-f33fdc9b552f", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 0:\n", "3 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1Bin 2
0A10.00.0100.00.0
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1 Bin 2\n", "0 A 10.0 0.0 100.0 0.0" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio = BioSim1D(n_bins=3, chem_data=chem_data)\n", "bio.inject_conc_to_bin(bin_address=1, delta_conc=100., chem_index=0)\n", "bio.describe_state()" ] }, { "cell_type": "code", "execution_count": 11, "id": "0fd086b3-8b8f-48cf-958f-50799c130bd3", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 0.033333:\n", "3 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1Bin 2
0A10.033.33333.33433.333
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1 Bin 2\n", "0 A 10.0 33.333 33.334 33.333" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio.time_step_threshold = 0.34 # TO slightly BYPASS THE SAFETY LIMIT \n", " # (for the allowed maximum delta Time) that is typically in place\n", "\n", "bio.diffuse(time_step=0.033333, n_steps=1)\n", "bio.describe_state()" ] }, { "cell_type": "markdown", "id": "df3c48e7-7492-499c-a4ac-d3a2f1e9dcaa", "metadata": {}, "source": [ "#### When the time step is (0.33333 / diffusion rate),\n", "#### a 3-bin system, configured as above, equilibrates in *a single step!*\n", "#### Again, an excessively large time step! (note that we're using the default value of 1 for delta_x)" ] }, { "cell_type": "markdown", "id": "25d75d27-675a-458b-a465-4b94b610a627", "metadata": {}, "source": [ "## Start over again with a new 3-bin system : same as the last round, but this time use a delta_x = 10 (rather than the default 1)" ] }, { "cell_type": "code", "execution_count": 12, "id": "8659eb72-6304-4512-b452-6acd84cc23af", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 0:\n", "3 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1Bin 2
0A10.00.0100.00.0
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1 Bin 2\n", "0 A 10.0 0.0 100.0 0.0" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio = BioSim1D(n_bins=3, chem_data=chem_data)\n", "bio.inject_conc_to_bin(bin_address=1, delta_conc=100., chem_index=0)\n", "bio.describe_state()" ] }, { "cell_type": "code", "execution_count": 13, "id": "e55a0df6-0442-44d1-a749-cbb356c53e58", "metadata": { "lines_to_next_cell": 2, "tags": [] }, "outputs": [], "source": [ "bio.time_step_threshold = 0.34 # TO slightly BYPASS THE SAFETY LIMIT \n", " # (for the allowed maximum delta Time) that is typically in place\n", " \n", "# When the time step is (delta_x**2 * 0.33333 / diffusion rate),\n", "# a 3-bin system, configured as above, equilibrates in a single step!" ] }, { "cell_type": "code", "execution_count": 14, "id": "a93fc736-b8a7-4442-9ff5-cf9d622e7259", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "SYSTEM STATE at Time t = 3.3333:\n", "3 bins and 1 chemical species\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SpeciesDiff rateBin 0Bin 1Bin 2
0A10.033.33333.33433.333
\n", "
" ], "text/plain": [ " Species Diff rate Bin 0 Bin 1 Bin 2\n", "0 A 10.0 33.333 33.334 33.333" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bio.diffuse(time_step=3.3333, n_steps=1, delta_x=10)\n", "bio.describe_state()" ] }, { "cell_type": "markdown", "id": "0d369711-0001-41b4-bbf0-f0b8964317d1", "metadata": {}, "source": [ "### This generalizes the previous result to any arbitrary delta_x :\n", "#### in this scenario, a time step of (delta_x**2 * 0.33333 / diffusion rate)\n", "#### lead to equilibrium in a single step!" ] }, { "cell_type": "code", "execution_count": null, "id": "64fc8359-10c2-46f9-affa-6a9e048c22b4", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "78d68b19-4031-421a-b030-186ede007bef", "metadata": {}, "source": [ "# The physical intuition explored in this experiment suggests enforcing:\n", "#### `delta_t < delta_x**2 * 0.33333 / diffusion rate`\n", "\n", "#### Interestingly, this is only slightly stricter than the upper bound that emerges from the **\"Von Neumann stability analysis\"** of the diffusion equation in 1-D, which states that solutions may become unstable unless\n", "delta_t < delta_x**2 * 0.5 / diffusion rate\n", "\n", "## NOTE: The method `BioSim1D.diffuse_step()` enforces the stricter upper bound" ] }, { "cell_type": "markdown", "id": "91ff8434-597e-4fd6-b76f-7f3572732863", "metadata": {}, "source": [ "#### To keep in mind that the _\"Von Neumann stability analysis\"_ is ONLY applicable when the diffusion equation is being solved with the **\"explicit Forward-Time Centered Space\"** method, which is the one currently used by the `diffuse_step()` method. \n", "An explanation can be found at: https://www.youtube.com/watch?v=QUiUGNwNNmo" ] }, { "cell_type": "code", "execution_count": null, "id": "08fb1f3a", "metadata": {}, "outputs": [], "source": [] } ], "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.13" } }, "nbformat": 4, "nbformat_minor": 5 }