Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Andrew cheng827 01 vib (WIP) #96

Open
wants to merge 7 commits into
base: master
Choose a base branch
from
Open

Andrew cheng827 01 vib (WIP) #96

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
33b3f70
Corrected \mathcal{I} to I. This is a test commit for Andrew
AndrewCheng827 Jul 12, 2023
34ff50a
Fixed typos in chapter 1; Transformed 'solver' functions using Devito…
AndrewCheng827 Jul 19, 2023
f6936f8
Chapter 1 now passes nbval tests
AndrewCheng827 Aug 1, 2023
644fc08
Passed all pytests on chapters 2,3, and 4
AndrewCheng827 Aug 7, 2023
ec3df8b
Chapter 1 (WIP)
AndrewCheng827 Sep 7, 2023
919d874
WIP
AndrewCheng827 Sep 14, 2023
8474b23
WIP
AndrewCheng827 Sep 14, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
141 changes: 86 additions & 55 deletions fdm-devito-notebooks/01_vib/vib_app.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,29 +2,19 @@
"cells": [
{
"cell_type": "code",
"execution_count": 3,
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Notebook initialized with ipy backend.\n"
"ename": "ModuleNotFoundError",
"evalue": "No module named 'mayavi'",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mModuleNotFoundError\u001b[0m Traceback (most recent call last)",
"Cell \u001b[0;32mIn[1], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mmayavi\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m mlab\n\u001b[1;32m 2\u001b[0m mlab\u001b[38;5;241m.\u001b[39minit_notebook()\n\u001b[1;32m 3\u001b[0m mlab\u001b[38;5;241m.\u001b[39mtest_plot3d()\n",
"\u001b[0;31mModuleNotFoundError\u001b[0m: No module named 'mayavi'"
]
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "de8153dff7fe4161829af04aaeff07bf",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Image(value=b'\\x89PNG\\r\\n\\x1a\\n\\x00\\x00\\x00\\rIHDR\\x00\\x00\\x01\\x90\\x00\\x00\\x01^\\x08\\x02\\x00\\x00\\x00$?\\xde_\\x00\\…"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
Expand Down Expand Up @@ -55,7 +45,7 @@
"\n",
"The following text derives some of the most well-known physical\n",
"problems that lead to second-order ODE models of the type addressed in\n",
"this ${DOCUMENT}. We consider a simple spring-mass system; thereafter\n",
"this book. We consider a simple spring-mass system; thereafter\n",
"extended with nonlinear spring, damping, and external excitation; a\n",
"spring-mass system with sliding friction; a simple and a physical\n",
"(classical) pendulum; and an elastic pendulum.\n",
Expand All @@ -75,22 +65,22 @@
"\n",
"The most fundamental mechanical vibration system is depicted in [Figure](#vib:app:mass_spring:fig). A body with mass $m$ is attached to a\n",
"spring and can move horizontally without friction (in the wheels). The\n",
"position of the body is given by the vector $\\rpos(t) = u(t)\\ii$, where\n",
"$\\ii$ is a unit vector in $x$ direction.\n",
"position of the body is given by the vector $\\textbf{r}(t) = u(t)\\textbf{i}$, where\n",
"$\\textbf{i}$ is a unit vector in $x$ direction.\n",
"There is\n",
"only one force acting on the body: a spring force $\\F_s =-ku\\ii$, where\n",
"only one force acting on the body: a spring force $\\textbf{F}_s =-ku\\textbf{i}$, where\n",
"$k$ is a constant. The point $x=0$, where $u=0$, must therefore\n",
"correspond to the body's position\n",
"where the spring is neither extended nor compressed, so the force\n",
"vanishes.\n",
"\n",
"The basic physical principle that governs the motion of the body is\n",
"Newton's second law of motion: $\\F=m\\acc$, where\n",
"$\\F$ is the sum of forces on the body, $m$ is its mass, and $\\acc=\\ddot\\rpos$\n",
"Newton's second law of motion: $\\textbf{F}=m\\textbf{a}$, where\n",
"$\\textbf{F}$ is the sum of forces on the body, $m$ is its mass, and $\\textbf{a}=\\ddot{r}$\n",
"is the acceleration. We use the dot for differentiation with respect\n",
"to time, which is\n",
"usual in mechanics. Newton's second law simplifies here\n",
"to $-\\F_s=m\\ddot u\\ii$, which translates to"
"to $-\\textbf{F}_s=m\\ddot u\\textbf{i}$, which translates to"
]
},
{
Expand Down Expand Up @@ -129,7 +119,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"mathcal{I}_t is\n",
"It is\n",
"not uncommon to divide by $m$\n",
"and introduce the frequency $\\omega = \\sqrt{k/m}$:"
]
Expand All @@ -156,7 +146,7 @@
"This is the model problem in the first part of this chapter, with the\n",
"small difference that we write the time derivative of $u$ with a dot\n",
"above, while we used $u^{\\prime}$ and $u^{\\prime\\prime}$ in previous\n",
"parts of the ${DOCUMENT}.\n",
"parts of the book.\n",
"\n",
"\n",
"Since only one scalar mathematical quantity, $u(t)$, describes the\n",
Expand Down Expand Up @@ -239,18 +229,18 @@
"The mechanical system in [Figure](#vib:app:mass_spring:fig) can easily be\n",
"extended to the more general system in [Figure](#vib:app:mass_gen:fig),\n",
"where the body is attached to a spring and a dashpot, and also subject\n",
"to an environmental force $F(t)\\ii$. The system has still only one\n",
"to an environmental force $F(t)\\textbf{i}$. The system has still only one\n",
"degree of freedom since the body can only move back and forth parallel to\n",
"the $x$ axis. The spring force was linear, $\\F_s=-ku\\ii$,\n",
"the $x$ axis. The spring force was linear, $\\textbf{F}_s=-ku\\textbf{i}$,\n",
"in the section [Oscillating mass attached to a spring](#vib:app:mass_spring), but in more general cases it can\n",
"depend nonlinearly on the position. We therefore set $\\F_s=s(u)\\ii$.\n",
"depend nonlinearly on the position. We therefore set $\\textbf{F}_s=s(u)\\textbf{i}$.\n",
"The dashpot, which acts\n",
"as a damper, results in a force $\\F_d$ that depends on the body's\n",
"velocity $\\dot u$ and that always acts against the motion.\n",
"The mathematical model of the force is written $\\F_d =f(\\dot u)\\ii$.\n",
"A positive $\\dot u$ must result in a force acting in the positive $x$\n",
"as a damper, results in a force $\\textbf{F}_d$ that depends on the body's\n",
"velocity $\\dot{u}$ and that always acts against the motion.\n",
"The mathematical model of the force is written $\\textbf{F}_d =f(\\dot u)\\textbf{i}$.\n",
"A positive $\\dot{u}$ must result in a force acting in the positive $x$\n",
"direction.\n",
"Finally, we have the external environmental force $\\F_e = F(t)\\ii$.\n",
"Finally, we have the external environmental force $\\textbf{F}_e = F(t)\\textbf{i}$.\n",
"\n",
"Newton's second law of motion now involves three forces:"
]
Expand All @@ -260,7 +250,7 @@
"metadata": {},
"source": [
"$$\n",
"F(t)\\ii - f(\\dot u)\\ii - s(u)\\ii = m\\ddot u \\ii\\thinspace .\n",
"F(t)\\textbf{i} - f(\\dot u)\\textbf{i} - s(u)\\textbf{i} = m\\ddot{u} \\textbf{i} \\thinspace .\n",
"$$"
]
},
Expand Down Expand Up @@ -786,7 +776,6 @@
"motion of the pendulum *and* the size of the forces during the motion.\n",
"The present section exemplifies how to make such a dynamic body\n",
"diagram.\n",
"% if FORMAT == 'pdflatex':\n",
"Two typical snapshots of free body diagrams are displayed below\n",
"(the drag force is magnified 5 times to become more visual!).\n",
"\n",
Expand All @@ -798,17 +787,44 @@
"\n",
"<!-- end figure -->\n",
"\n",
"\n",
"% else:\n",
"<!-- dom:MOVIE: [https://github.com/hplgit/pysketcher/raw/master/doc/pub/tutorial/mov-tut/pendulum/movie.mp4] The drag force is magnified 5 times! % endif -->\n",
"<!-- begin movie -->"
]
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": 4,
"metadata": {},
"outputs": [],
"outputs": [
{
"data": {
"text/html": [
"\n",
"<div>\n",
"<video loop controls width='640' height='365' preload='none'>\n",
" <source src='https://github.com/hplgit/pysketcher/raw/master/doc/pub/tutorial/mov-tut/pendulum/movie.mp4' type='video/mp4; codecs=\"avc1.42E01E, mp4a.40.2\"'>\n",
" <source src='https://github.com/hplgit/pysketcher/raw/master/doc/pub/tutorial/mov-tut/pendulum/movie.webm' type='video/webm; codecs=\"vp8, vorbis\"'>\n",
" <source src='https://github.com/hplgit/pysketcher/raw/master/doc/pub/tutorial/mov-tut/pendulum/movie.ogg' type='video/ogg; codecs=\"theora, vorbis\"'>\n",
"</video>\n",
"</div>\n",
"<p><em>The drag force is magnified 5 times!</em></p>\n",
"\n",
"<!-- Issue warning if in a Safari browser -->\n",
"<script language=\"javascript\">\n",
"if (!!(window.safari)) {\n",
" document.write(\"<div style=\"width: 95%%; padding: 10px; border: 1px solid #100; border-radius: 4px;\"><p><font color=\"red\">The above movie will not play in Safari - use Chrome, Firefox, or Opera.</font></p></div>\")}\n",
"</script>\n",
"\n"
],
"text/plain": [
"<IPython.core.display.HTML object>"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from IPython.display import HTML\n",
"_s = \"\"\"\n",
Expand All @@ -819,7 +835,7 @@
" <source src='https://github.com/hplgit/pysketcher/raw/master/doc/pub/tutorial/mov-tut/pendulum/movie.ogg' type='video/ogg; codecs=\"theora, vorbis\"'>\n",
"</video>\n",
"</div>\n",
"<p><em>The drag force is magnified 5 times! % endif</em></p>\n",
"<p><em>The drag force is magnified 5 times!</em></p>\n",
"\n",
"<!-- Issue warning if in a Safari browser -->\n",
"<script language=\"javascript\">\n",
Expand Down Expand Up @@ -966,10 +982,12 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"! Should I convert this to Devito?\n",
"\n",
"def simulate(alpha, Theta, dt, T):\n",
" import odespy\n",
"\n",
Expand Down Expand Up @@ -1238,7 +1256,7 @@
"The position vector $\\rpos$ in Cartesian coordinates reads\n",
"$\\rpos(t) = x(t)\\ii + y(t)\\jj$, where $\\ii$ and $\\jj$ are unit vectors\n",
"in the $x$ and $y$ directions, respectively.\n",
"mathcal{I}_t is convenient to introduce the Cartesian components $n_x$ and $n_y$\n",
"It is convenient to introduce the Cartesian components $n_x$ and $n_y$\n",
"of the normal vector:"
]
},
Expand Down Expand Up @@ -2091,11 +2109,7 @@
"is found in a separate project on the web:\n",
"<https://github.com/hplgit/bumpy>. You may start looking at the\n",
"\"tutorial\":\n",
"% if FORMAT == 'pdflatex':\n",
"\"http://hplgit.github.io/bumpy/doc/pub/bumpy.pdf\".\n",
"% else:\n",
"\"http://hplgit.github.io/bumpy/doc/pub/bumpy.html\".\n",
"% endif\n",
"\n",
"## Bouncing ball\n",
"<div id=\"vib:app:bouncing_ball\"></div>\n",
Expand Down Expand Up @@ -2737,7 +2751,7 @@
"metadata": {},
"source": [
"$$\n",
"s(u) = \\frac{k}{\\alpha}\\tanh(\\alpha u) = ku + \\frac{1}{3}\\alpha^2 ku^3 + \\frac{2}{15}\\alpha^4 k u^5 + \\Oof{u^6}\\thinspace .\n",
"s(u) = \\frac{k}{\\alpha}\\tanh(\\alpha u) = ku + \\frac{1}{3}\\alpha^2 ku^3 + \\frac{2}{15}\\alpha^4 k u^5 + \\mathcal{O}({u^6})\\thinspace .\n",
"$$"
]
},
Expand Down Expand Up @@ -2929,6 +2943,8 @@
"metadata": {},
"outputs": [],
"source": [
"! Should I change this to Devito?\n",
"\n",
"def simulate(beta, gamma, delta=0,\n",
" num_periods=8, time_steps_per_period=60):\n",
" # Use oscillations without friction to set dt and T\n",
Expand Down Expand Up @@ -3029,10 +3045,21 @@
},
{
"cell_type": "code",
"execution_count": 11,
"execution_count": 6,
"metadata": {},
"outputs": [],
"outputs": [
{
"ename": "SyntaxError",
"evalue": "Missing parentheses in call to 'print'. Did you mean print(...)? (4188502967.py, line 42)",
"output_type": "error",
"traceback": [
"\u001b[0;36m Cell \u001b[0;32mIn[6], line 42\u001b[0;36m\u001b[0m\n\u001b[0;31m print '%4d v=%8.5f h=%8.5f %s' % (n, v[n+1], h[n+1], mode)\u001b[0m\n\u001b[0m ^\u001b[0m\n\u001b[0;31mSyntaxError\u001b[0m\u001b[0;31m:\u001b[0m Missing parentheses in call to 'print'. Did you mean print(...)?\n"
]
}
],
"source": [
"! Devito-ize this\n",
"\n",
"import numpy as np\n",
"\n",
"def solver(H, C_R, dt, T, eps_v=0.01, eps_h=0.01):\n",
Expand Down Expand Up @@ -3626,6 +3653,8 @@
"metadata": {},
"outputs": [],
"source": [
"! Devito-ize the odespy part?\n",
"\n",
"import odespy\n",
"import numpy as np\n",
"import scitools.std as plt\n",
Expand Down Expand Up @@ -4231,6 +4260,8 @@
"metadata": {},
"outputs": [],
"source": [
"! Devito-ize odespy\n",
"\n",
"def simulate_drag(\n",
" beta=0.9, # dimensionless elasticity parameter\n",
" gamma=0, # dimensionless drag parameter\n",
Expand Down Expand Up @@ -4369,7 +4400,7 @@
"metadata": {},
"source": [
"**a)**\n",
"mathcal{I}_t is easy to show that $x(t)$ and $y(t)$ go like sine and cosine\n",
"It is easy to show that $x(t)$ and $y(t)$ go like sine and cosine\n",
"functions. Use this idea to derive the exact solution.\n",
"\n",
"\n",
Expand Down Expand Up @@ -5433,7 +5464,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -5447,7 +5478,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.3"
"version": "3.10.6"
}
},
"nbformat": 4,
Expand Down
Loading
Loading