From a6129b0bdf96cee5011801bfe69c8ba06e96e85e Mon Sep 17 00:00:00 2001 From: kangwonlee Date: Tue, 28 May 2024 04:17:34 +0000 Subject: [PATCH] Clean ipynb 8bfa77e185580035b1d83788455e6e584d655f35 --- 50_ode/30_Runge_Kutta.ipynb | 1779 ++++++++++++++++------------------- 1 file changed, 821 insertions(+), 958 deletions(-) diff --git a/50_ode/30_Runge_Kutta.ipynb b/50_ode/30_Runge_Kutta.ipynb index 2abf68aa..1b32b94a 100644 --- a/50_ode/30_Runge_Kutta.ipynb +++ b/50_ode/30_Runge_Kutta.ipynb @@ -1,960 +1,823 @@ { - "cells": [ - { - "cell_type": "markdown", - "metadata": { - "id": "view-in-github", - "colab_type": "text" - }, - "source": [ - "\"Open" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "Y7kIXItMXbS6" - }, - "outputs": [], - "source": [ - "# This cell is for the Google Colaboratory\n", - "# https://stackoverflow.com/a/63519730\n", - "if 'google.colab' in str(get_ipython()):\n", - " path_py = '/content/nmisp_py'\n", - "\n", - " import os\n", - " if not os.path.exists(path_py):\n", - " import subprocess\n", - " subprocess.run(\n", - " ('git', 'clone', 'https://github.com/kangwonlee/nmisp_py')\n", - " )\n", - " assert os.path.exists(path_py)\n", - "\n", - " import sys\n", - " sys.path.insert(0, path_py)\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "A1lRFYItXbS-" - }, - "outputs": [], - "source": [ - "# 그래프, 수학 기능 추가\n", - "# Add graph and math features\n", - "import pylab as py\n", - "import numpy as np\n", - "import numpy.linalg as nl\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "8b8WbTLsXbS-" - }, - "source": [ - "# 룽게-쿠타법 (RK4)
Runge-Kutta Method (RK4)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "fuzBJPPGXbS_" - }, - "source": [ - "오일러법과 훈의 방법은 $t=t_0$, $t=t_1$, ... 에서의 기울기만 사용하였다.
Euler Method and Heun's Method used slopes at $t=t_0$, $t=t_1$, ... only.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "A5aHLE9WXbTA" - }, - "source": [ - "룽게-쿠타 법은 1900년대에 독일 수학자 칼 룽게와 마틴 쿠타가 공동 개발했던 상미분 방정식의 근사 해법의 모음이다.
Runge-Kutta methods are a group of numerical methods to solve ordinary differential equations that two German mathematicians Carl Runge and Martin Kutta developed in 1900s.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "EI2B4REKXbTA" - }, - "source": [ - "실은 오일러법이나 훈의 방법도 룽게-쿠타법에 포함된다.
In fact, Euler method or Heun's Method are included in the Runge-Kutta method.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "9lVTakrHXbTB" - }, - "source": [ - "룽게-쿠타법 가운데 **RK4**가 가장 널리 사용된다.
Among the Runge-Kutta methods, **RK4** is used the most frequently.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "lC7azkhjXbTC" - }, - "source": [ - "**RK4** 는 $t=t_1$ 에서의 해를 구하기 위해 $t=t_0$, $t=t_\\frac{1}{2}=\\frac{1}{2}\\left(t_0+t_1\\right)$, $t=t_1$ 에서의 기울기를 사용한다.
\n", - "To find a solution at $t=t_1$, **RK4** uses slopes at $t=t_0$, $t=t_\\frac{1}{2}=\\frac{1}{2}\\left(t_0+t_1\\right)$, and $t=t_1$.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "rpWz44ofXbTD" - }, - "source": [ - "$$\n", - "\\begin{cases}\n", - " \\frac{d}{dt}\\mathbf{x}=f(t, \\mathbf{x}) \\\\\n", - " \\mathbf{x}(t_0)=\\mathbf{x}_0\n", - "\\end{cases}\n", - "$$\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "6il3HHA0XbTE" - }, - "source": [ - "위 미분방정식의 경우, 시간 간격 $\\Delta t=t_1 - t_0$ 일 때 절차는 다음과 같다.
\n", - "For the differential equation above, when time step $\\Delta t=t_1 - t_0$, the procedure is as follows.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "NBDuJNFWXbTE" - }, - "source": [ - "* RK4\n", - "\n", - "| Equation 수식 | Description 설명 |\n", - "|:-------------------:|:-----------------------------------------------:|\n", - "| $$\\mathbf{s}_1=f(t_0, \\mathbf{x}_0)$$ | $t_0$ 에서의 기울기 $\\mathbf{s}_1$을 구한다
At $t_0$, find the slope $\\mathbf{s}_1$ |\n", - "| $$\\mathbf{s}_2=f(t_{\\frac{1}{2}}, \\mathbf{x}_0 + \\mathbf{s}_1 {\\frac{1}{2}} \\Delta t)$$ | $t_0$ 로 부터 ${\\frac{1}{2}} \\Delta t$ 만큼 기울기 $\\mathbf{s}_1$을 따라 전진하여 $t_\\frac{1}{2}$ 에서의 기울기 $\\mathbf{s}_2$을 구한다.
Advancing from $t_0$ by ${\\frac{1}{2}} \\Delta t$ in time, along the slope $\\mathbf{s}_1$, find the slope$\\mathbf{s}_2$ at $t_\\frac{1}{2}$ |\n", - "| $$ \\mathbf{s}_3=f(t_{\\frac{1}{2}}, \\mathbf{x}_0 + \\mathbf{s}_2 {\\frac{1}{2}} \\Delta t) $$ | 다시 한번, $t_0$에서 $t_{\\frac{1}{2}}$ 까지 $\\mathbf{s}_2$를 따라 전진하여 $t_\\frac{1}{2}$ 에서의 기울기 $\\mathbf{s}_3$을 구한다.
Once again, by time-advancing from $t_0$ to $t=t_\\frac{1}{2}$, along the slope $\\mathbf{s}_2$, find the slope $\\mathbf{s}_3$. |\n", - "| $$ \\mathbf{s}_4=f(t_1, \\mathbf{x}_0 + \\mathbf{s}_3 \\Delta t) $$ |이번에는 $t_0$에서 $t_1$ 까지 $\\mathbf{s}_3$를 따라 전진하여 $t=t_1$ 에서의 기울기 $\\mathbf{s}_4$을 구한다.
This time, by going forward from $t_0$ to $t_1$, find the slope at $t=t_1$, $\\mathbf{s}_4$ |\n", - "| $$ \\mathbf{s}_{RK4}=\\frac{1}{6} \\left( \\mathbf{s}_1 + 2\\mathbf{s}_2 + 2\\mathbf{s}_3 + \\mathbf{s}_4 \\right) $$ | $t_0 \\le t \\le t_1$ 구간을 대표하는 기울기 $\\mathbf{s}_{RK4}$을 구한다.
Find the slope $\\mathbf{s}_{RK4}$ representing interval $t_0 \\le t \\le t_1$.|\n", - "| $$ \\mathbf{x}_1 = \\mathbf{x}_0 + \\mathbf{s}_{RK4} \\Delta t $$ | $t_0$ 로 부터 $t_1$ 까지 $\\mathbf{s}_{RK4}$를 따라 전진하여 $\\mathbf{x}_1$ 을 구한다.
Find $\\mathbf{x}_1$ by advacning from $t_0$ to $t_1$ following slope $\\mathbf{s}_{RK4}$. |\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "FIy7PzBYXbTE" - }, - "source": [ - "| Euler | Heun | RK4 |\n", - "|:------:|:------:|:------:|\n", - "| $$ s_1 = f(t_0, x_0) $$| $$ s_1 = f(t_0, x_0) $$ | $$s_1=f(t_0, x_0)$$ |\n", - "| $$ x_1=x_0 + s_1 \\Delta t $$ | $$ s_2=f(t_{i+1}, x_0 + s_1 \\Delta t) t $$ | $$ s_2=f(t_{\\frac{1}{2}}, x_0 + s_1 {\\frac{1}{2}} \\Delta t) $$ |\n", - "| $$ $$ | $$ s_{Heun} = \\frac{1}{2}(s_1 + s_2) $$ | $$ s_3=f(t_{\\frac{1}{2}}, x_0 + s_2 {\\frac{1}{2}} \\Delta t) $$ |\n", - "| $$ $$ | $$ x_1 = x_0 + s_{Heun} \\Delta t $$ | $$ s_4=f(t_1, x_0 + s_3 \\Delta t) $$ |\n", - "| | | $$ s_{RK4}=\\frac{1}{6} \\left( s_1 + 2s_2 + 2s_3 + s_4 \\right) $$ |\n", - "| | | $$ x_1 = x_0 + s_{RK4} \\Delta t $$ |\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "hczdVTcFXbTF" - }, - "source": [ - "python 으로 써 보자.
Let's write in python.\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "yGKh8YlEXbTF" - }, - "outputs": [], - "source": [ - "def rk4_step(f, x0, t0, t1):\n", - " \"\"\"\n", - " One time step of Runge-Kutta method\n", - "\n", - " f: dx_dt function\n", - " x0 : initial condition\n", - " t0 : this step time\n", - " t1 : next step time\n", - " \"\"\"\n", - " delta_t = (t1 - t0)\n", - " delta_t_half = delta_t * 0.5\n", - " t_half = t0 + delta_t_half\n", - "\n", - " # Step 1\n", - " s1 = f(t0, x0)\n", - "\n", - " # Step 2\n", - " s2 = f(t_half, x0 + s1 * delta_t_half)\n", - "\n", - " # Step 3\n", - " s3 = f(t_half, x0 + s2 * delta_t_half)\n", - "\n", - " # Step 4\n", - " s4 = f(t1, x0 + s3 * delta_t)\n", - "\n", - " # Step 5\n", - " s = (1.0 / 6.0) * (s1 + (s2 + s3) * 2 + s4)\n", - "\n", - " # Step 6\n", - " x1 = x0 + s * delta_t\n", - "\n", - " return x1\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "QHbndoBCXbTF" - }, - "outputs": [], - "source": [ - "def rk4(f, t_array, x_0):\n", - " time_list = [t_array[0]]\n", - " result_list = [x_0]\n", - "\n", - " x_i = x_0\n", - "\n", - " for k, t_i in enumerate(t_array[:-1]):\n", - " # time step\n", - " x_i_plus_1 = rk4_step(f, x_i, t_i, t_array[k+1])\n", - "\n", - " time_list.append(t_array[k+1])\n", - " result_list.append(x_i_plus_1)\n", - "\n", - " x_i = x_i_plus_1\n", - "\n", - " return time_list, result_list\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "KyPoPVawXbTG" - }, - "source": [ - "다시 1계 선형 상미분 방정식의 예를 살펴 보자.
Let's reconsider the 1st order linear ODE.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "t0SnaDMLXbTG" - }, - "source": [ - "$$\n", - "\\left\\{\n", - " \\begin{align}\n", - " a_0 \\frac{d}{dt}x(t)+a_1 x(t)&=0 \\\\\n", - " x(0)&=x_0 \\\\\n", - " \\end{align}\n", - "\\right.\n", - "$$\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "6Jynb0b4XbTG" - }, - "source": [ - "룽게-쿠타법 결과를 엄밀해, 오일러법, 훈법과 비교해보자.
Let's compare the result from Runge-Kutta method with the exact solution, Euler method, and Heun's method.\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "uxhU12laXbTG" - }, - "outputs": [], - "source": [ - "a_0, a_1 = 2.0, 1.0\n", - "\n", - "def dx_dt(t, x):\n", - " return - a_1 * x / a_0\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "ON9vvYWtXbTH" - }, - "source": [ - "Initial value $x(t_0)$
초기값 $x(t_0)$\n", - "\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "Igp-ki4EXbTH" - }, - "outputs": [], - "source": [ - "x_0 = 4.5\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "-IUZr6HfXbTG" - }, - "outputs": [], - "source": [ - "def exact(t):\n", - " return x_0 * py.exp((-a_1 / a_0) * t)\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "SJnzP_h6XbTH" - }, - "outputs": [], - "source": [ - "import ode_solver\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "SHcvu03sXbTH" - }, - "source": [ - "$\\Delta t$\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "MXKYrKWtXbTH" - }, - "outputs": [], - "source": [ - "delta_t = 0.01\n", - "t_sec_array = np.arange(0, 6 + delta_t*0.5, delta_t)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "nu26ptz_XbTH" - }, - "source": [ - "오일러법
Euler method\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "sM8c9bf0XbTH" - }, - "outputs": [], - "source": [ - "t_euler_out, x_euler_out = ode_solver.euler(dx_dt, t_sec_array, x_0)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "4vMxm9NhXbTI" - }, - "source": [ - "훈의 방법
\n", - "Heun's method\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "S00_GdXEXbTI" - }, - "outputs": [], - "source": [ - "t_heun__out, x_heun__out = ode_solver.heun(dx_dt, t_sec_array, x_0)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "C4VsspuUXbTI" - }, - "source": [ - "RK4\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "2Cue7ZsGXbTI" - }, - "outputs": [], - "source": [ - "t_rk4___out, x_rk4___out = rk4(dx_dt, t_sec_array, x_0)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "4TqIuv-bXbTI" - }, - "source": [ - "근사해 그림
\n", - "Plots of approximate solutions\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "c0UwCF3KXbTI" - }, - "outputs": [], - "source": [ - "import ode_plot\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "scrolled": false, - "id": "1fSUSji4XbTI" - }, - "outputs": [], - "source": [ - "py.figure(figsize=(16, 9))\n", - "\n", - "# Approximate solutions\n", - "py.plot(t_euler_out, x_euler_out, '.-', label='Euler')\n", - "py.plot(t_heun__out, x_heun__out, '.-', label='Heun')\n", - "py.plot(t_rk4___out, x_rk4___out, '.-', label='RK4')\n", - "\n", - "# *** Exact Solution\n", - "t_exact_array = np.linspace(0, 6)\n", - "exact = ode_plot.ExactPlotterFirstOrderODE(\n", - " t_exact_array,\n", - " a_0=a_0, a_1=a_1, x_0=x_0)\n", - "\n", - "exact.plot()\n", - "\n", - "py.xlabel('t(sec)')\n", - "py.ylabel('x(m)')\n", - "\n", - "py.legend(loc=0)\n", - "py.grid(True)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "mHmpMjYhXbTJ" - }, - "source": [ - "룽게-쿠타법의 해가 엄밀해에 더 가까움을 알 수 있다.
\n", - "We can see that Runge-Kutta method is closer to the exact solution.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "SBbez0ozXbTJ" - }, - "source": [ - "## Scipy\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "j21PYNZUXbTJ" - }, - "outputs": [], - "source": [ - "import scipy.integrate as si\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "bgDS0IueXbTJ" - }, - "outputs": [], - "source": [ - "sol = si.solve_ivp(dx_dt, (t_heun__out[0], t_heun__out[-1]), [x_0], t_eval=t_heun__out)\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "9YxBzrKMXbTJ" - }, - "outputs": [], - "source": [ - "py.figure(figsize=(16, 9))\n", - "\n", - "py.plot(sol.t, sol.y[0, :], 'o', label='solve_ivp')\n", - "py.plot(t_euler_out, x_euler_out, '.-', label='Euler')\n", - "py.plot(t_heun__out, x_heun__out, '*-', label='Heun')\n", - "py.plot(t_rk4___out, x_rk4___out, '.-', label='RK4')\n", - "\n", - "# plot exact solution\n", - "exact = ode_plot.ExactPlotterFirstOrderODE(\n", - " py.array(t_rk4___out),\n", - " a_0=a_0, a_1=a_1, x_0=x_0)\n", - "exact.plot()\n", - "\n", - "py.grid(True)\n", - "py.xlabel('t(sec)')\n", - "py.ylabel('y(t)')\n", - "py.legend(loc=0);\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "vdwSqHJ5XbTN" - }, - "outputs": [], - "source": [ - "import pandas as pd\n", - "\n", - "\n", - "df = pd.DataFrame(\n", - " data={\n", - " 'euler':x_euler_out,\n", - " 'heun' :x_heun__out,\n", - " 'rk4' :x_rk4___out,\n", - " 'solve_ivp':sol.y[0, :],\n", - " 'exact':exact.exact(py.array(t_heun__out))\n", - " },\n", - " index=pd.Series(t_heun__out, name='t(sec)'),\n", - " columns=['exact', 'euler', 'heun', 'rk4', 'solve_ivp']\n", - ")\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "XJWHXwV4XbTN" - }, - "outputs": [], - "source": [ - "df['euler_error'] = df.euler - df.exact\n", - "df['heun_error'] = df.heun - df.exact\n", - "df['rk4_error'] = df.rk4 - df.exact\n", - "df['solve_ivp_error'] = df.solve_ivp - df.exact\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "dAwtl9LkXbTN" - }, - "source": [ - "표 형태
Table form\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "ZiJ15iwXXbTN" - }, - "outputs": [], - "source": [ - "pd.set_option('display.max_rows', 10)\n", - "df\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "NwPWQ1HBXbTN" - }, - "source": [ - "각종 통계
Statistics\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "zayLZbZuXbTN" - }, - "outputs": [], - "source": [ - "df.describe()\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "ghkHxDLkXbTN" - }, - "source": [ - "이 경우, RK4 오차에 대한 의견은?
\n", - "In this case, what do you think about the error of the RK4?\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "ABaSIH1kXbTN" - }, - "outputs": [], - "source": [ - "import numpy.linalg as nl\n", - "\n", - "\n", - "nl.norm(df.euler_error), nl.norm(df.heun_error), nl.norm(df.rk4_error), nl.norm(df.solve_ivp_error),\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "kwhGp05dXbTO" - }, - "source": [ - "### 계산시간
Computation time\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "SzsWXiyjXbTO" - }, - "source": [ - "$\\Delta t$ 을 더 작은 값으로 바꾸어 보자.
Let's try a smaller $\\Delta t$.\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "31_p5vQNXbTO" - }, - "outputs": [], - "source": [ - "delta_t = 1e-3\n", - "t_sec_array = np.arange(0, 6 + delta_t*0.5, delta_t)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "mEFL9xWpXbTO" - }, - "source": [ - "오일러법
Euler method\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "VWD-W4RlXbTO" - }, - "outputs": [], - "source": [ - "%%timeit -n100\n", - "\n", - "t_euler_out, x_euler_out = ode_solver.euler(dx_dt, t_sec_array, x_0)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "Ud6aM3kzXbTO" - }, - "source": [ - "훈의 방법
\n", - "Heun's method\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "Zi8sulUPXbTO" - }, - "outputs": [], - "source": [ - "%%timeit -n100\n", - "\n", - "t_heun__out, x_heun__out = ode_solver.heun(dx_dt, t_sec_array, x_0)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "InuHpF-yXbTO" - }, - "source": [ - "RK4\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "IQTOYiHZXbTP" - }, - "outputs": [], - "source": [ - "%%timeit -n100\n", - "\n", - "t_rk4___out, x_rk4___out = rk4(dx_dt, t_sec_array, x_0)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "XrOMpkSsXbTP" - }, - "source": [ - "계산의 단계 수가 더 많은 해법일 수록 계산 시간이 많이 필요한 것으로 보인다.
\n", - "With more steps, computation takes more time.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "va0CrDbTXbTP" - }, - "source": [ - "## 연습 문제
Exercises\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "k5BqschNXbTP" - }, - "source": [ - "다음 미분방정식의 엄밀해를 구하시오:
\n", - "Find exact solution of the following differential equation:\n", - "\n", - "$$\n", - "\\begin{align}\n", - "10 \\frac{d}{dt}x(t) + 50 x(t) &= 0 \\\\\n", - "x(0) &= 2\n", - "\\end{align}\n", - "$$\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "wzP_tQ1WXbTP" - }, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "m2QCDZr3XbTP" - }, - "source": [ - "위 미분방정식의 수치해를 오일러법으로 구하시오.
\n", - "Find numerical solution of the above differential equation using Euler Method.\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "LS0uFQZBXbTP" - }, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "I1yAKjSXXbTP" - }, - "source": [ - "위 미분방정식의 수치해를 훈의 방법으로 구하고 엄밀해, 오일러법과 비교하시오.
\n", - "Find numerical solution of the above differential equation using Heun's method and compare with exact solution and Euler Method.\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "t78hZToGXbTQ" - }, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "nH3Vi54KXbTQ" - }, - "source": [ - "위 미분방정식의 수치해를 RK4법으로 구하고 엄밀해, 오일러법, 훈의 방법과 비교하시오.
\n", - "Find numerical solution of the above differential equation using RK$ and compare with exact solution, Euler Method, and Heun's Method.\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "1gzqvCswXbTQ" - }, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "ePJNR10cXbTQ" - }, - "source": [ - "## Final Bell
마지막 종\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "v4yV1_TTXbTQ" - }, - "outputs": [], - "source": [ - "# stackoverfow.com/a/24634221\n", - "import os\n", - "os.system(\"printf '\\a'\");\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "id": "GfRJP_CgXbTQ" - }, - "outputs": [], - "source": [] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Python 3", - "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.7.3" - }, - "colab": { - "provenance": [], - "include_colab_link": true - } - }, - "nbformat": 4, - "nbformat_minor": 0 + "cells": [ + { + "cell_type": "markdown", + "metadata": { + "id": "view-in-github", + "colab_type": "text" + }, + "source": [ + "\"Open\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# This cell is for the Google Colaboratory\n", + "# https://stackoverflow.com/a/63519730\n", + "if 'google.colab' in str(get_ipython()):\n", + " path_py = '/content/nmisp_py'\n", + "\n", + " import os\n", + " if not os.path.exists(path_py):\n", + " import subprocess\n", + " subprocess.run(\n", + " ('git', 'clone', 'https://github.com/kangwonlee/nmisp_py')\n", + " )\n", + " assert os.path.exists(path_py)\n", + "\n", + " import sys\n", + " sys.path.insert(0, path_py)\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# 그래프, 수학 기능 추가\n", + "# Add graph and math features\n", + "import pylab as py\n", + "import numpy as np\n", + "import numpy.linalg as nl\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# 룽게-쿠타법 (RK4)
Runge-Kutta Method (RK4)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "오일러법과 훈의 방법은 $t=t_0$, $t=t_1$, ... 에서의 기울기만 사용하였다.
Euler Method and Heun's Method used slopes at $t=t_0$, $t=t_1$, ... only.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "룽게-쿠타 법은 1900년대에 독일 수학자 칼 룽게와 마틴 쿠타가 공동 개발했던 상미분 방정식의 근사 해법의 모음이다.
Runge-Kutta methods are a group of numerical methods to solve ordinary differential equations that two German mathematicians Carl Runge and Martin Kutta developed in 1900s.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "실은 오일러법이나 훈의 방법도 룽게-쿠타법에 포함된다.
In fact, Euler method or Heun's Method are included in the Runge-Kutta method.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "룽게-쿠타법 가운데 **RK4**가 가장 널리 사용된다.
Among the Runge-Kutta methods, **RK4** is used the most frequently.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**RK4** 는 $t=t_1$ 에서의 해를 구하기 위해 $t=t_0$, $t=t_\\frac{1}{2}=\\frac{1}{2}\\left(t_0+t_1\\right)$, $t=t_1$ 에서의 기울기를 사용한다.
\n", + "To find a solution at $t=t_1$, **RK4** uses slopes at $t=t_0$, $t=t_\\frac{1}{2}=\\frac{1}{2}\\left(t_0+t_1\\right)$, and $t=t_1$.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "$$\n", + "\\begin{cases}\n", + " \\frac{d}{dt}\\mathbf{x}=f(t, \\mathbf{x}) \\\\\n", + " \\mathbf{x}(t_0)=\\mathbf{x}_0\n", + "\\end{cases}\n", + "$$\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "위 미분방정식의 경우, 시간 간격 $\\Delta t=t_1 - t_0$ 일 때 절차는 다음과 같다.
\n", + "For the differential equation above, when time step $\\Delta t=t_1 - t_0$, the procedure is as follows.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "* RK4\n", + "\n", + "| Equation 수식 | Description 설명 |\n", + "|:-------------------:|:-----------------------------------------------:|\n", + "| $$\\mathbf{s}_1=f(t_0, \\mathbf{x}_0)$$ | $t_0$ 에서의 기울기 $\\mathbf{s}_1$을 구한다
At $t_0$, find the slope $\\mathbf{s}_1$ |\n", + "| $$\\mathbf{s}_2=f(t_{\\frac{1}{2}}, \\mathbf{x}_0 + \\mathbf{s}_1 {\\frac{1}{2}} \\Delta t)$$ | $t_0$ 로 부터 ${\\frac{1}{2}} \\Delta t$ 만큼 기울기 $\\mathbf{s}_1$을 따라 전진하여 $t_\\frac{1}{2}$ 에서의 기울기 $\\mathbf{s}_2$을 구한다.
Advancing from $t_0$ by ${\\frac{1}{2}} \\Delta t$ in time, along the slope $\\mathbf{s}_1$, find the slope$\\mathbf{s}_2$ at $t_\\frac{1}{2}$ |\n", + "| $$ \\mathbf{s}_3=f(t_{\\frac{1}{2}}, \\mathbf{x}_0 + \\mathbf{s}_2 {\\frac{1}{2}} \\Delta t) $$ | 다시 한번, $t_0$에서 $t_{\\frac{1}{2}}$ 까지 $\\mathbf{s}_2$를 따라 전진하여 $t_\\frac{1}{2}$ 에서의 기울기 $\\mathbf{s}_3$을 구한다.
Once again, by time-advancing from $t_0$ to $t=t_\\frac{1}{2}$, along the slope $\\mathbf{s}_2$, find the slope $\\mathbf{s}_3$. |\n", + "| $$ \\mathbf{s}_4=f(t_1, \\mathbf{x}_0 + \\mathbf{s}_3 \\Delta t) $$ |이번에는 $t_0$에서 $t_1$ 까지 $\\mathbf{s}_3$를 따라 전진하여 $t=t_1$ 에서의 기울기 $\\mathbf{s}_4$을 구한다.
This time, by going forward from $t_0$ to $t_1$, find the slope at $t=t_1$, $\\mathbf{s}_4$ |\n", + "| $$ \\mathbf{s}_{RK4}=\\frac{1}{6} \\left( \\mathbf{s}_1 + 2\\mathbf{s}_2 + 2\\mathbf{s}_3 + \\mathbf{s}_4 \\right) $$ | $t_0 \\le t \\le t_1$ 구간을 대표하는 기울기 $\\mathbf{s}_{RK4}$을 구한다.
Find the slope $\\mathbf{s}_{RK4}$ representing interval $t_0 \\le t \\le t_1$.|\n", + "| $$ \\mathbf{x}_1 = \\mathbf{x}_0 + \\mathbf{s}_{RK4} \\Delta t $$ | $t_0$ 로 부터 $t_1$ 까지 $\\mathbf{s}_{RK4}$를 따라 전진하여 $\\mathbf{x}_1$ 을 구한다.
Find $\\mathbf{x}_1$ by advacning from $t_0$ to $t_1$ following slope $\\mathbf{s}_{RK4}$. |\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "| Euler | Heun | RK4 |\n", + "|:------:|:------:|:------:|\n", + "| $$ s_1 = f(t_0, x_0) $$| $$ s_1 = f(t_0, x_0) $$ | $$s_1=f(t_0, x_0)$$ |\n", + "| $$ x_1=x_0 + s_1 \\Delta t $$ | $$ s_2=f(t_{i+1}, x_0 + s_1 \\Delta t) t $$ | $$ s_2=f(t_{\\frac{1}{2}}, x_0 + s_1 {\\frac{1}{2}} \\Delta t) $$ |\n", + "| $$ $$ | $$ s_{Heun} = \\frac{1}{2}(s_1 + s_2) $$ | $$ s_3=f(t_{\\frac{1}{2}}, x_0 + s_2 {\\frac{1}{2}} \\Delta t) $$ |\n", + "| $$ $$ | $$ x_1 = x_0 + s_{Heun} \\Delta t $$ | $$ s_4=f(t_1, x_0 + s_3 \\Delta t) $$ |\n", + "| | | $$ s_{RK4}=\\frac{1}{6} \\left( s_1 + 2s_2 + 2s_3 + s_4 \\right) $$ |\n", + "| | | $$ x_1 = x_0 + s_{RK4} \\Delta t $$ |\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "python 으로 써 보자.
Let's write in python.\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def rk4_step(f, x0, t0, t1):\n", + " \"\"\"\n", + " One time step of Runge-Kutta method\n", + "\n", + " f: dx_dt function\n", + " x0 : initial condition\n", + " t0 : this step time\n", + " t1 : next step time\n", + " \"\"\"\n", + " delta_t = (t1 - t0)\n", + " delta_t_half = delta_t * 0.5\n", + " t_half = t0 + delta_t_half\n", + "\n", + " # Step 1\n", + " s1 = f(t0, x0)\n", + "\n", + " # Step 2\n", + " s2 = f(t_half, x0 + s1 * delta_t_half)\n", + "\n", + " # Step 3\n", + " s3 = f(t_half, x0 + s2 * delta_t_half)\n", + "\n", + " # Step 4\n", + " s4 = f(t1, x0 + s3 * delta_t)\n", + "\n", + " # Step 5\n", + " s = (1.0 / 6.0) * (s1 + (s2 + s3) * 2 + s4)\n", + "\n", + " # Step 6\n", + " x1 = x0 + s * delta_t\n", + "\n", + " return x1\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def rk4(f, t_array, x_0):\n", + " time_list = [t_array[0]]\n", + " result_list = [x_0]\n", + "\n", + " x_i = x_0\n", + "\n", + " for k, t_i in enumerate(t_array[:-1]):\n", + " # time step\n", + " x_i_plus_1 = rk4_step(f, x_i, t_i, t_array[k+1])\n", + "\n", + " time_list.append(t_array[k+1])\n", + " result_list.append(x_i_plus_1)\n", + "\n", + " x_i = x_i_plus_1\n", + "\n", + " return time_list, result_list\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "다시 1계 선형 상미분 방정식의 예를 살펴 보자.
Let's reconsider the 1st order linear ODE.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "$$\n", + "\\left\\{\n", + " \\begin{align}\n", + " a_0 \\frac{d}{dt}x(t)+a_1 x(t)&=0 \\\\\n", + " x(0)&=x_0 \\\\\n", + " \\end{align}\n", + "\\right.\n", + "$$\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "룽게-쿠타법 결과를 엄밀해, 오일러법, 훈법과 비교해보자.
Let's compare the result from Runge-Kutta method with the exact solution, Euler method, and Heun's method.\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "a_0, a_1 = 2.0, 1.0\n", + "\n", + "def dx_dt(t, x):\n", + " return - a_1 * x / a_0\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Initial value $x(t_0)$
초기값 $x(t_0)$\n", + "\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "x_0 = 4.5\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def exact(t):\n", + " return x_0 * py.exp((-a_1 / a_0) * t)\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import ode_solver\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "$\\Delta t$\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "delta_t = 0.01\n", + "t_sec_array = np.arange(0, 6 + delta_t*0.5, delta_t)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "오일러법
Euler method\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "t_euler_out, x_euler_out = ode_solver.euler(dx_dt, t_sec_array, x_0)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "훈의 방법
\n", + "Heun's method\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "t_heun__out, x_heun__out = ode_solver.heun(dx_dt, t_sec_array, x_0)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "RK4\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "t_rk4___out, x_rk4___out = rk4(dx_dt, t_sec_array, x_0)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "근사해 그림
\n", + "Plots of approximate solutions\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import ode_plot\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "scrolled": false + }, + "outputs": [], + "source": [ + "py.figure(figsize=(16, 9))\n", + "\n", + "# Approximate solutions\n", + "py.plot(t_euler_out, x_euler_out, '.-', label='Euler')\n", + "py.plot(t_heun__out, x_heun__out, '.-', label='Heun')\n", + "py.plot(t_rk4___out, x_rk4___out, '.-', label='RK4')\n", + "\n", + "# *** Exact Solution\n", + "t_exact_array = np.linspace(0, 6)\n", + "exact = ode_plot.ExactPlotterFirstOrderODE(\n", + " t_exact_array,\n", + " a_0=a_0, a_1=a_1, x_0=x_0)\n", + "\n", + "exact.plot()\n", + "\n", + "py.xlabel('t(sec)')\n", + "py.ylabel('x(m)')\n", + "\n", + "py.legend(loc=0)\n", + "py.grid(True)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "룽게-쿠타법의 해가 엄밀해에 더 가까움을 알 수 있다.
\n", + "We can see that Runge-Kutta method is closer to the exact solution.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Scipy\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import scipy.integrate as si\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "sol = si.solve_ivp(dx_dt, (t_heun__out[0], t_heun__out[-1]), [x_0], t_eval=t_heun__out)\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "py.figure(figsize=(16, 9))\n", + "\n", + "py.plot(sol.t, sol.y[0, :], 'o', label='solve_ivp')\n", + "py.plot(t_euler_out, x_euler_out, '.-', label='Euler')\n", + "py.plot(t_heun__out, x_heun__out, '*-', label='Heun')\n", + "py.plot(t_rk4___out, x_rk4___out, '.-', label='RK4')\n", + "\n", + "# plot exact solution\n", + "exact = ode_plot.ExactPlotterFirstOrderODE(\n", + " py.array(t_rk4___out),\n", + " a_0=a_0, a_1=a_1, x_0=x_0)\n", + "exact.plot()\n", + "\n", + "py.grid(True)\n", + "py.xlabel('t(sec)')\n", + "py.ylabel('y(t)')\n", + "py.legend(loc=0);\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd\n", + "\n", + "\n", + "df = pd.DataFrame(\n", + " data={\n", + " 'euler':x_euler_out,\n", + " 'heun' :x_heun__out,\n", + " 'rk4' :x_rk4___out,\n", + " 'solve_ivp':sol.y[0, :],\n", + " 'exact':exact.exact(py.array(t_heun__out))\n", + " },\n", + " index=pd.Series(t_heun__out, name='t(sec)'),\n", + " columns=['exact', 'euler', 'heun', 'rk4', 'solve_ivp']\n", + ")\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "df['euler_error'] = df.euler - df.exact\n", + "df['heun_error'] = df.heun - df.exact\n", + "df['rk4_error'] = df.rk4 - df.exact\n", + "df['solve_ivp_error'] = df.solve_ivp - df.exact\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "표 형태
Table form\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "pd.set_option('display.max_rows', 10)\n", + "df\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "각종 통계
Statistics\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "df.describe()\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "이 경우, RK4 오차에 대한 의견은?
\n", + "In this case, what do you think about the error of the RK4?\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy.linalg as nl\n", + "\n", + "\n", + "nl.norm(df.euler_error), nl.norm(df.heun_error), nl.norm(df.rk4_error), nl.norm(df.solve_ivp_error),\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 계산시간
Computation time\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "$\\Delta t$ 을 더 작은 값으로 바꾸어 보자.
Let's try a smaller $\\Delta t$.\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "delta_t = 1e-3\n", + "t_sec_array = np.arange(0, 6 + delta_t*0.5, delta_t)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "오일러법
Euler method\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "%%timeit -n100\n", + "\n", + "t_euler_out, x_euler_out = ode_solver.euler(dx_dt, t_sec_array, x_0)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "훈의 방법
\n", + "Heun's method\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "%%timeit -n100\n", + "\n", + "t_heun__out, x_heun__out = ode_solver.heun(dx_dt, t_sec_array, x_0)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "RK4\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "%%timeit -n100\n", + "\n", + "t_rk4___out, x_rk4___out = rk4(dx_dt, t_sec_array, x_0)\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "계산의 단계 수가 더 많은 해법일 수록 계산 시간이 많이 필요한 것으로 보인다.
\n", + "With more steps, computation takes more time.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 연습 문제
Exercises\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "다음 미분방정식의 엄밀해를 구하시오:
\n", + "Find exact solution of the following differential equation:\n", + "\n", + "$$\n", + "\\begin{align}\n", + "10 \\frac{d}{dt}x(t) + 50 x(t) &= 0 \\\\\n", + "x(0) &= 2\n", + "\\end{align}\n", + "$$\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "위 미분방정식의 수치해를 오일러법으로 구하시오.
\n", + "Find numerical solution of the above differential equation using Euler Method.\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "위 미분방정식의 수치해를 훈의 방법으로 구하고 엄밀해, 오일러법과 비교하시오.
\n", + "Find numerical solution of the above differential equation using Heun's method and compare with exact solution and Euler Method.\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "위 미분방정식의 수치해를 RK4법으로 구하고 엄밀해, 오일러법, 훈의 방법과 비교하시오.
\n", + "Find numerical solution of the above differential equation using RK$ and compare with exact solution, Euler Method, and Heun's Method.\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Final Bell
마지막 종\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# stackoverfow.com/a/24634221\n", + "import os\n", + "os.system(\"printf '\\a'\");\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "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.7.3" + }, + "colab": { + "provenance": [], + "include_colab_link": true + } + }, + "nbformat": 4, + "nbformat_minor": 0 } \ No newline at end of file