Skip to content

Commit

Permalink
update to solutions for tutorial 7
Browse files Browse the repository at this point in the history
  • Loading branch information
mjksill committed Nov 1, 2023
1 parent f093744 commit 171b1d6
Showing 1 changed file with 115 additions and 2 deletions.
117 changes: 115 additions & 2 deletions tutorials/tutorialS1_7-soln.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -203,13 +203,126 @@
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "a4a6f792-8aee-4208-9c9d-aa7a14a13d76",
"metadata": {},
"source": [
"To determine the mean concentration $\\bar{c}$ at the end of the pipe, we first need to determine the total amount of the solute in a small disk of length $\\Delta z$ at the end of the pipe. The total amount of solute in the pipe is\n",
"$$\n",
"\\begin{align*}\n",
"\\Delta N &= \\int_0^R \\Delta z 2\\pi r dr\\, c(r)\n",
"\\end{align*}\n",
"$$\n",
"To get the mean concentration, we need to divide the total amount of solute in the disk by the volume of the disk, which is $\\Delta V = \\Delta z \\pi R^2$. This gives\n",
"$$\n",
"\\begin{align*}\n",
"\\bar{c} &= \\frac{\\Delta N}{\\Delta V}\n",
"\\\\\n",
"&= \\frac{1}{\\Delta z \\pi R^2} \\int_0^R \\Delta z 2\\pi r dr\\, c(r)\n",
"\\\\\n",
"&= \\frac{2}{R^2} \\int_0^R dr\\, r c(r)\n",
"\\end{align*}\n",
"$$\n",
"\n",
"\n",
"To evaluate this integral, we can use the trapezoid rule code we developed in Q1 of tutorial 1 or the code from Q2 of tutorial 5. "
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d4ddc5f1",
"metadata": {},
"outputs": [],
"source": []
"source": [
"import numpy as np\n",
"import pylab as plt\n",
"\n",
"def v(r):\n",
" v0 = 1.0e-2\n",
" x = r/R\n",
" return 2.0*v0*(1.0-x**2)\n",
"\n",
"def c(r):\n",
" c0 = 0.2\n",
" L = 1.0\n",
" t = k*L/v(r)\n",
" #print(t)\n",
" return c0*np.exp(-t)\n",
"\n",
"\n",
"R = 1.0e-3\n",
"\n",
"\n",
"### Below is the code from Q2 of tutorial 5\n",
"### First we need to modify the function we are integrating:\n",
"def I_f(r):\n",
" return r*c(r)\n",
"\n",
"N = 100+1\n",
"r_data = np.linspace(0, 0.9999*R, N) # this needs to change to match our limits of integration\n",
" # note that there is a problem exactly at r=R, so we will just\n",
" # integrate a little bit below.\n",
"\n",
"I = 0.0\n",
"for i in range(N-1):\n",
" ra = r_data[i]\n",
" rb = r_data[i+1]\n",
" I += 0.5*(rb-ra)*(I_f(ra)+I_f(rb))\n",
"\n",
"\n",
"print(I)\n",
"c_bar = 2/R**2 * I\n",
"print(f'c_bar = {c_bar:.5f} M')\n",
"\n"
]
},
{
"cell_type": "markdown",
"id": "a19612a1-4850-4662-aee2-bd172ac5b368",
"metadata": {},
"source": [
"We don't need to write the trapezoid rule code for ourselves. There are many routines available on the internet that can do numerical integration for us. One of them is the `quad` function from the `scipy` library."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "907e3901-b66c-4862-8aed-9881dffc1c3d",
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import pylab as plt\n",
"\n",
"def v(r):\n",
" v0 = 1.0e-2\n",
" x = r/R\n",
" return 2.0*v0*(1.0-x**2)\n",
"\n",
"def c(r):\n",
" c0 = 0.2\n",
" L = 1.0\n",
" t = k*L/v(r)\n",
" #print(t)\n",
" return c0*np.exp(-t)\n",
"\n",
"\n",
"\n",
"R = 1.0e-3\n",
"\n",
"\n",
"from scipy.integrate import quad\n",
"I = quad(I_f, 0, R) # the output of quad is a tuple. The first element is the value of the integral. The\n",
" # second number is the estimated uncertainty of the value.\n",
"print(I)\n",
"c_bar = 2/R**2 * I[0] # note: we use I[0], because we just want to use sthe value of the integral.\n",
"print(f'c_bar = {c_bar:.5f} M')\n",
"\n",
"\n",
"\n"
]
},
{
"cell_type": "markdown",
Expand Down Expand Up @@ -1289,7 +1402,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.8"
"version": "3.11.5"
}
},
"nbformat": 4,
Expand Down

0 comments on commit 171b1d6

Please sign in to comment.