Skip to content
Projects
Groups
Snippets
Help
Loading...
Sign in / Register
Toggle navigation
N
numer-60-2
Project
Project
Details
Activity
Cycle Analytics
Repository
Repository
Files
Commits
Branches
Tags
Contributors
Graph
Compare
Charts
Issues
0
Issues
0
List
Board
Labels
Milestones
Merge Requests
0
Merge Requests
0
CI / CD
CI / CD
Pipelines
Jobs
Schedules
Charts
Wiki
Wiki
Snippets
Snippets
Members
Members
Collapse sidebar
Close sidebar
Activity
Graph
Charts
Create a new issue
Jobs
Commits
Issue Boards
Open sidebar
Littichai Buddaken
numer-60-2
Commits
a67cb5bd
Commit
a67cb5bd
authored
Apr 05, 2018
by
Littichai Buddaken
Browse files
Options
Browse Files
Download
Email Patches
Plain Diff
add file Numerical-Differentiation-2
parent
6b8215e5
Hide whitespace changes
Inline
Side-by-side
Showing
1 changed file
with
509 additions
and
0 deletions
+509
-0
Numerical-Differentiation-2.ipynb
week09/Numerical-Differentiation-2.ipynb
+509
-0
No files found.
week09/Numerical-Differentiation-2.ipynb
0 → 100644
View file @
a67cb5bd
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Numerical Differentiation\n",
"การหาผลต่าง \n",
"\n",
"> กำหนดฟังก์ชัน $f(x)$ จงคำนวณหา $\\frac{d^nf}{dx^n}$ เมื่อกำหนดค่า $x$ ให้\n",
"\n",
"โดยที่ $f(x)$ กำหนดได้สองแบบคือ\n",
"\n",
"* กำหนดเป็นฟังก์ชันสำหรับ การคำนวณหาโดยตรง $y = f(x)$ \n",
"* กำหนดความสัมพันธ์เป็นแบบชุดข้อมูลมาให้ $(x_i,y_i), i=0,1,..,n$\n",
"\n",
"## Taylor Series of $f(x)$ at $a$ is\n",
"$$\n",
"\\sum_{n=0}^{\\infty}\\frac{f^n(a)}{n!}(x-a)^n\n",
"$$\n",
"\n",
"ตัวอย่าง Taylor series สำหรับ $e^x$ ที่ $a = 0$ คือ\n",
"$$\n",
"\\frac{x^0}{0!} + \\frac{x^1}{1!} + \\frac{x^2}{2!} + \\frac{x^3}{3!} + ... \\\\\n",
"1 + x + \\frac{x^2}{2!} + \\frac{x^3}{3!} + ... \\\\\n",
"\\sum_{n=0}^{\\infty}\\frac{x^n}{n!}\n",
"$$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Backward/Forward Taylor Series on $f(x)$\n",
"$$\n",
"f(x+h) = f(x) + hf^1(x) + \\frac{h^2}{2!}f^2(x) + \\frac{h^3}{3!}f^3(x) + \\frac{h^4}{4!}f^4(x) + ... \\to (a) \\\\\n",
"f(x-h) = f(x) - hf^1(x) + \\frac{h^2}{2!}f^2(x) - \\frac{h^3}{3!}f^3(x) + \\frac{h^4}{4!}f^4(x) - ... \\to (b) \\\\\n",
"f(x+2h) = f(x) + 2hf^1(x) + \\frac{(2h)^2}{2!}f^2(x) + \\frac{(2h)^3}{3!}f^3(x) + \\frac{(2h)^4}{4!}f^4(x) + ... \\to (c) \\\\\n",
"f(x-2h) = f(x) - 2hf^1(x) + \\frac{(2h)^2}{2!}f^2(x) - \\frac{(2h)^3}{3!}f^3(x) + \\frac{(2h)^4}{4!}f^4(x) - ... \\to (d) \\\\\n",
"$$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### ผลต่างของสมการ $(a), (b), (c), (d)$\n",
"$ (a) + (b) $\n",
"$$\n",
"f(x+h)+f(x-h) = 2f(x) + h^2f^2(x) + \\frac{h^4}{12}f^4(x) + ... \\to (e) \\\\\n",
"$$\n",
"$ (a) - (b) $\n",
"$$\n",
"f(x+h)-f(x-h) = 2hf^1(x) + \\frac{h^3}{3}f^3(x) + ... \\to (f) \\\\\n",
"$$\n",
"$ (c) + (d) $\n",
"$$\n",
"f(x+2h)+f(x-2h) = 2f(x) + 4h^2f^2(x) + \\frac{4h^4}{3}f^4(x) + ... \\to (g) \\\\\n",
"$$\n",
"$ (c) - (d) $\n",
"$$\n",
"f(x+2h)-f(x-2h) = 4hf^1(x) + \\frac{8h^3}{3}f^3(x) + ... \\to (h) \\\\\n",
"$$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Finite Difference Approximations\n",
"\n",
"### Central Difference Approximations: CDA\n",
"\n",
"#### First CDA\n",
"$$\n",
"f^1(x) = \\frac{f(x+h) - f(x-h)}{2h} - \\frac{h^2}{3!}f^3(x) - \\frac{h^4}{5!}f^5(x) - \\frac{h^6}{7!}f^7(x) -...\n",
"$$\n",
"หรือ\n",
"$$\n",
"f'(x) = \\frac{f(x+h)-f(x-h)}{2h} + O(h^2)\n",
"$$\n",
"สูตร\n",
"$$\n",
"f'(x) \\approx \\frac{f(x+h)-f(x-h)}{2h} \\to (c)\n",
"$$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Exercise\n",
"1. จงคำนวณหาค่าประมาณของ $f^1(x)$, (First CDA) พร้อมค่าผิดพลาด $O(h^2)$ ของฟังก์ชันต่อไปนี้ เมื่อกำหนด $h = 0.001, x=1$\n",
" * $f(x) = e^{-x}$\n",
" * $f(x) = sin(x)e^{-cos(x)}$\n",
"\n",
"2. จงเขียนฟังก์ชันเพื่อแสดงตารางค่าความสัมพันธ์ระหว่าง $h$ กับ ค่า First CDA ของฟังก์ชันต่างๆในข้อ 1. เมื่อกำหนด $x = 1$\n",
" และ $h = [0.64, 0.32, 0.16, 0.08, 0.04, 0.02, 0.01, 0.005, 0.0025, 0.00125]$\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Second CDA\n",
"$$\n",
"f^2(x) = \\frac{f(x+h) - 2f(x) + f(x-h)}{h^2} + O(h^2)\n",
"$$\n",
"สูตร\n",
"$$\n",
"f^2(x) \\approx \\frac{f(x+h) - 2f(x) + f(x-h)}{h^2} \\to (d)\n",
"$$\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Exercise\n",
"1. จงคำนวณหาค่าประมาณของ $f^2(x)$, (Second CDA) พร้อมค่าผิดพลาด $O(h^2)$ ของฟังก์ชันต่อไปนี้ เมื่อกำหนด $h = 0.001, x=1$\n",
" * $f(x) = e^{-x}$\n",
" * $f(x) = sin(x)e^{-cos(x)}$\n",
"\n",
"2. จงเขียนฟังก์ชันเพื่อแสดงตารางค่าความสัมพันธ์ระหว่าง $h$ กับ ค่า Second CDA ของฟังก์ชันต่างๆในข้อ 1. เมื่อกำหนด $x = 1$\n",
" และ $h = [0.64, 0.32, 0.16, 0.08, 0.04, 0.02, 0.01, 0.005, 0.0025, 0.00125]$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Third CDA\n",
"$$\n",
"f^3(x) = \\frac{f(x+2h) - 2f(x+h) + 2f(x-h) + -f(x-2h)}{2h^3} + O(h^2) \n",
"$$\n",
"สูตร\n",
"$$\n",
"f^3(x) \\approx \\frac{f(x+2h) - 2f(x+h) + 2f(x-h) + -f(x-2h)}{2h^3} \\to (e)\n",
"$$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Exercise\n",
"1. จงคำนวณหาค่าประมาณของ $f^1(x)$, (Third CDA) พร้อมค่าผิดพลาด $O(h^2)$ ของฟังก์ชันต่อไปนี้ เมื่อกำหนด $h = 0.001, x=1$\n",
" * $f(x) = e^{-x}$\n",
" * $f(x) = sin(x)e^{-cos(x)}$\n",
"\n",
"2. จงเขียนฟังก์ชันเพื่อแสดงตารางค่าความสัมพันธ์ระหว่าง $h$ กับ ค่า Third CDA ของฟังก์ชันต่างๆในข้อ 1. เมื่อกำหนด $x = 1$\n",
" และ $h = [0.64, 0.32, 0.16, 0.08, 0.04, 0.02, 0.01, 0.005, 0.0025, 0.00125]$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Fourth CDA\n",
"$$\n",
"f^4(x) = \\frac{f(x+2h) - 4f(x+h) + 6f(x) - 4f(x-h) + f(x-2h)}{h^4} + O(h^2)\n",
"$$\n",
"สูตร\n",
"$$\n",
"f^4(x) \\approx \\frac{f(x+2h) - 4f(x+h) + 6f(x) - 4f(x-h) + f(x-2h)}{h^4} \\to (f)\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Exercise\n",
"1. จงคำนวณหาค่าประมาณของ $f^1(x)$, (Fourth CDA) พร้อมค่าผิดพลาด $O(h^2)$ ของฟังก์ชันต่อไปนี้ เมื่อกำหนด $h = 0.001, x=1$\n",
" * $f(x) = e^{-x}$\n",
" * $f(x) = sin(x)e^{-cos(x)}$\n",
"\n",
"2. จงเขียนฟังก์ชันเพื่อแสดงตารางค่าความสัมพันธ์ระหว่าง $h$ กับ ค่า Fourth CDA ของฟังก์ชันต่างๆในข้อ 1. เมื่อกำหนด $x = 1$\n",
" และ $h = [0.64, 0.32, 0.16, 0.08, 0.04, 0.02, 0.01, 0.005, 0.0025, 0.00125]$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Example\n",
"\n",
"1. กำหนดฟังก์ชัน $f(x) = e^{-x}$ จงคำนวณหา $f^2(x)$ ที่ $x = 1, h = 0.01$ โดยใช้ second CDA\n",
"> จากสูตร $(d)$ เราสามารถประมาณค่า $f^2(x)$ ได้ โดยใช้สมการ\n",
"$$\n",
"f^2(x) = \\frac{f(x-h) - 2f(x) + f(x+h)}{h^2}\n",
"$$\n",
"ดังนั้น คำตอบคือ\n",
"$$\n",
"f^2(1) = \\frac{f(1-0.01) - 2\\times f(1) + f(1+0.01)}{0.01^2} \\\\\n",
"f^2(1) = \\frac{f(0.98) - 2\\times f(1) + f(1.01)}{0.0001} \\\\\n",
"f^2(1) = \\frac{e^{-0.98} - 2\\times e^{-1} + e^{-1.01}}{0.0001} \\\\\n",
"$$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Exercise\n",
"\n",
"1. กำหนดฟังก์ชัน $f(x) = x^3e^x$ จงคำนวณหา $f^2(x)$ ที่ $x = 1, h = 0.01$ โดยใช้ second CDA\n",
"\n",
"2. กำหนดฟังก์ชัน $f(x) = x^2e^{-x}$ จงคำนวณหา $f^3(x)$ ที่ $x = 1, h = 0.01$ โดยใช้ third CDA\n",
"\n",
"3. กำหนดฟังก์ชัน $f(x) = \\frac{e^{-x}}{x^3}$ จงคำนวณหา $f^4(x)$ ที่ $x = 1, h = 0.01$ โดยใช้ fourth CDA\n",
"\n",
"4. กำหนดฟังก์ชัน $f(x) = \\frac{e^{-x}}{x^3}$ จงเขียนโปรแกรมเพื่อคำนวณหา $f^4(x)$ ที่ $x = 1, h = 0.01$ โดยใช้ fourth CDA"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# สรุป Central Difference Approximations\n",
"\n",
"$$\n",
"(1) \\to f'(x) = \\frac{f(x+h)-f(x-h)}{2h} + O(h^2) \\\\\n",
"(2) \\to f^2(x) = \\frac{f(x+h) - 2f(x) + f(x-h)}{h^2} + O(h^2) \\\\\n",
"(3) \\to f^3(x) = \\frac{f(x+2h) - 2f(x+h) + 2f(x-h) + -f(x-2h)}{2h^3} + O(h^2) \\\\\n",
"(4) \\to f^4(x) = \\frac{f(x+2h) - 4f(x+h) + 6f(x) - 4f(x-h) + f(x-2h)}{h^4} + O(h^2) \\\\\n",
"$$\n",
"\n",
"หรือถ้าไม่คิด $O(h^2)$\n",
"\n",
"$(1) \\to f'(x) = \\frac{f(x+h)-f(x-h)}{2h}$\n",
"\n",
"$(2) \\to f^2(x) = \\frac{f(x+h) - 2f(x) + f(x-h)}{h^2}$\n",
"\n",
"$(3) \\to f^3(x) = \\frac{f(x+2h) - 2f(x+h) + 2f(x-h) + -f(x-2h)}{2h^3}$\n",
"\n",
"$(4) \\to f^4(x) = \\frac{f(x+2h) - 4f(x+h) + 6f(x) - 4f(x-h) + f(x-2h)}{h^4}$\n",
"\n",
"ปรับใหม่ \n",
"\n",
"$(1) \\to {2h}f'(x) = f(x+h)-f(x-h)$\n",
"\n",
"$(2) \\to {h^2}f^2(x) = f(x+h) - 2f(x) + f(x-h)$\n",
"\n",
"$(3) \\to {2h^3}f^3(x) = f(x+2h) - 2f(x+h) + 2f(x-h) + -f(x-2h)$\n",
"\n",
"$(4) \\to {h^4}f^4(x) = f(x+2h) - 4f(x+h) + 6f(x) - 4f(x-h) + f(x-2h)$\n",
"\n",
"แล้วเขียนออกมาเป็น ตาราง coefficient ได้เป็น\n",
"\n",
"# ตารางค่าสัมประสิทธิของ CDA $O(h^2$\n",
"\n",
"| สมการด้านซ้ายมือ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ |\n",
"|---------------------|----------------|---------------|-----------------|----------------|-----------------|\n",
"| | $f(x-2h)$ | $f(x-h)$ | $f(x)$ | $f(x+h)$ | $f(x+2h)$ |\n",
"|$2hf^1(x)$ | | -1 | 0 | 1 | |\n",
"|$h^2f^2(x)$ | | 1 | -2 | 1 | |\n",
"|$2h^3f^3(x)$ | -1 | 2 | 0 | -2 | 1 |\n",
"|$h^4f^4(x)$ | 1 | -4 | 6 | -4 | 1 |"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"ในกรณีที่ ให้มาเป็นจุด หรือไม่สามารถประเมินโดยใช้ค่า $+h, -h$ ได้ จะใช้\n",
"\n",
"\n",
"# ตารางค่าสัมประสิทธิของ Forward Finite Difference Approximation (FDA) $O(h^2$\n",
"\n",
"| สมการด้านซ้ายมือ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ |\n",
"|----------------|----------------|----------------|---------------|-----------------|----------------|-----------------|\n",
"| | $f(x)$ | $f(x+h)$ | $f(x+2h)$ | $f(x+3h)$ | $f(x+4h)$ | $f(x+5h)$ |\n",
"|$2hf^1(x)$ | -3 | 4 | -1 | | | |\n",
"|$h^2f^2(x)$ | 2 | -5 | 4 | -1 | | |\n",
"|$2h^3f^3(x)$ | -5 | 18 | -24 | 14 | -3 | |\n",
"|$h^4f^4(x)$ | 3 | -14 | 26 | -24 | 11 | -2 |\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"# ตารางค่าสัมประสิทธิของ Backward Finite Difference Approximation (FDA) $O(h^2$\n",
"\n",
"| สมการด้านซ้ายมือ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ | ค่าสัมประสิทธิ์ของ |\n",
"|----------------|----------------|----------------|---------------|-----------------|----------------|-----------------|\n",
"| | $f(x-5h)$ | $f(x-4h)$ | $f(x-3h)$ | $f(x-2h)$ | $f(x-h)$ | $f(x)$ |\n",
"|$2hf^1(x)$ | | | | 1 | -4 | 3 |\n",
"|$h^2f^2(x)$ | | | -1 | 4 | -5 | 2 |\n",
"|$2h^3f^3(x)$ | | 3 | -14 | 24 | 18 | 5 |\n",
"|$h^4f^4(x)$ | -2 | 11 | -24 | 26 | -14 | 3 |\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Richardson Extrapolation\n",
"\n",
"พิจารณา \n",
"$$\n",
"(2) \\to f^2(x) = \\frac{f(x+h) - 2f(x) + f(x-h)}{h^2} + O(h^2)\n",
"$$\n",
"ของฟังก์ชัน\n",
"$f(x) = e^{-x}$ ที่ $x=1$ \n",
"\n",
"เมื่อเปรียบเทียบกับค่าจริง(6 หลัก) $f^2(1) = e^{-1} = 0.3678794$\n",
"\n",
"\n",
"| $h$ | second CDA | ค่าจริง | $O(h^2)$ |\n",
"|---------|------------|----------|-----------|\n",
"| 0.64 | 0.380610 | 0.367879 | 0.012731 |\n",
"| 0.32 | 0.371035 | 0.367879 | 0.003156 |\n",
"| 0.16 | 0.368711 | 0.367879 | 0.000832 |\n",
"| 0.08 | 0.368281 | 0.367879 | 0.000402 |\n",
"| 0.04 | 0.36875 | 0.367879 | 0.000871 |\n",
"| 0.02 | 0.37 | 0.367879 | 0.002121 |\n",
"| 0.01 | 0.38 | 0.367879 | 0.012121 |\n",
"| 0.005 | 0.40 | 0.367879 | 0.032121 |\n",
"| 0.0025 | 0.48 | 0.367879 | 0.112121 |\n",
"| 0.00125 | 1.28 | 0.367879 | 0.912121 |\n"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"e**-1 = 0.36787944117144233\n"
]
}
],
"source": [
"import numpy as np\n",
"import math as ma\n",
"print('e**-1 =',np.e**-1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"เทคนิค Richardson Extrapolation ช่วยลดค่าความผิดพลาดได้ โดยมีเงื่อนไขคือ\n",
"$$\n",
"G = g(h) + E(h)\n",
"$$\n",
"และ $E(h) = ch^p$ เมื่อ $c,p$ เป็นค่าคงที่\n",
"\n",
"แทนค่าในสมการที่ $h_1, h_2$ \n",
"$$\n",
"G = g(h_1) + ch_{1}^p \\to (1) \\\\\n",
"G = g(h_2) + ch_{2}^p \\to (2)\n",
"$$\n",
"แล้วกำจัด $c$ จะได้\n",
"$$\n",
"G = \\frac{(h_1/h_2)^p g(h_2) - g(h_1)}{(h_1/h_2)^p -1}\n",
"$$\n",
"ซึ่งโดยปกติแล้วจะใช้ $h_2 = \\frac{h_1}{2}$ \n",
"จะทำให้\n",
"$$\n",
"G = \\frac{2^p g(h_1/2) - g(h_1)}{2^p - 1}\n",
"$$\n",
"\n",
"\n",
"## ตัวอย่าง\n",
"จากตาราง \n",
"$g(0.64) = 0.380610$ \n",
"$g(0.32) = 0.371035$\n",
"จะได้ $h_1 = 0.64, p=2$ \n",
"$$\n",
"G = \\frac{2^p g(h_1/2) - g(h_1)}{2^p - 1} \\\\\n",
"G = \\frac{2^2 g(0.64/2) - g(0.64)}{2^2 - 1} \\\\\n",
"G = \\frac{2^2 g(0.32) - g(0.64)}{2^2 - 1} \\\\\n",
"G = 0.367843\n",
"$$\n"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.38898750000000004"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"g= 2**2*0.371035-0.380610/2**2-1\n",
"g"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"def Extrapolation(g1,g2):\n",
" return ((2**2)*g2-g1)/2**2-1"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"-0.7241175"
]
},
"execution_count": 26,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Extrapolation(0.380610,0.371035)"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.36784333333333336"
]
},
"execution_count": 27,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def G(h):\n",
" cda={\n",
" 0.64:0.380610 ,\n",
" 0.32:0.371035 ,\n",
" 0.16:0.368711 ,\n",
" 0.08:0.368281 ,\n",
" 0.04:0.36875 ,\n",
" }\n",
" return cda[h]\n",
"def RichExtra(h,p):\n",
" return ((2**p)*G(h/2)-G(h))/(2**p -1)\n",
"RichExtra(0.64,2)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"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.6.1"
}
},
"nbformat": 4,
"nbformat_minor": 1
}
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment