มีเด็กๆที่สนใจคณิตศาสตร์มาถามผมว่าสมการทั่วๆไปแก้ออกมาเป๊ะๆไม่ได้แล้วเราทำอย่างไร ผมก็บอกว่าสมการส่วนใหญ่เราต้องหาคำตอบด้วยการประมาณเอาครับ ซึ่งวิธีอันหนึ่งที่เราสามารถใช้ได้ง่ายๆก็คือวิธีของนิวตัน ข้างล่างนี้เป็น
กระทู้ที่ผมเขียนไว้ที่
Mahidol Physics Educational Center ครับ:
******
สำหรับปัญหาที่เราแก้สมการโดยตรงไม่ได้ เราต้องแก้ด้วยวิธีประมาณด้วยตัวเลขครับ วิธีที่ใช้กันบ่อยๆวิธีหนึ่งก็คือวิธีการของนิวตัน (Newton's method:
http://en.wikipedia.org/wiki/Newton's_method) ครับ
วิธีการของนิวตันบอกว่า ถ้าจะแก้สมการ f(x) = 0 ให้เราเดาค่า x มาสักค่า (เรียกมันว่า x
0) ก็แล้วกัน แล้วค่า x อันต่อไป (เรียกมันว่า x
1) ที่น่าจะทำให้ f(x) ใกล้ศูนย์มากขึ้น ควรจะคำนวณอย่างนี้ครับ:
x
1 = x
0 - f(x
0)/f'(x
0) โดยที่ f'(x) คือ derivative ของ f(x) ครับ
ถ้าค่า x
1 ทำให้ f(x) ไม่ใกล้ 0 พอ เราก็หา x
2, x
3, x
4, ... ไปเรื่อยๆจนเราพอใจว่าค่า f(x
n) ใกล้ 0 พอแล้ว โดยที่ x
n หาได้จาก x
n-1 ดังนี้ครับ:
x
n = x
n-1 - f(x
n-1)/f'(x
n-1)
ถ้าจะทำการคำนวณด้วยวิธีของนิวตันใน
Mathematica สามารถทำอย่างนี้ครับ:
newtonSolve[f_, guess_, steps_] := NestList[ #1 - f[#1]/f'[#1] &, guess, steps]
f คือฟังค์ชั่นที่เราจะหา f(x) = 0
guess คือค่าที่เราเดาตอนแรกว่า f(guess) น่าจะไม่ห่างจาก 0 นัก
steps คือจำนวนครั้งที่เราจะทำการทำวิธีของนิวตันซำ้ๆกัน
ยกตัวอย่างเช่น เราจะหาค่ารูทที่สองของสอง เราก็เขียนสมการ f(x) = x^2 -2 = 0 ก่อน เพราะค่า x เท่ากับรูทที่สองของสองจะแก้สมการนั้นพอดี:
f[x_] := x^2 - 2
แล้วเราก็เรียก newtonSolve ด้วยฟังค์ชั่น f โดยเดาค่า guess = 1 และให้ทำซ้ำสักห้าครั้ง:
newtonSolve[f, 1, 5]
แล้วเราก็ได้ผลดังนี้: {1, 3/2, 17/12, 577/408, 665857/470832, 886731088897/627013566048}
Mathematica ทำการคำนวณให้เป็นค่าเศษส่วนไม่มีทศนิยม เพราะเราเดาด้วยค่า 1 ซึ่งเป็นจำนวนที่ไม่มีการประมาณเข้ามาเกี่ยวข้อง ถ้าเราต้องการคำตอบเป็นเลขทศนิยม เราก็ควรเดาด้วยค่า 1.0 ดังนี้:
newtonSolve[f, 1.0, 5]
แล้วเราก็จะได้ผลดังนี้: {1., 1.5, 1.41667, 1.41422, 1.41421, 1.41421} ซึ่งเราจะเห็นว่า 1.41421 นั้นเป็นค่าประมาณของรูทที่สองของสองได้ดีทีเดียว
นักศึกษาลองไปทดลองดูครับ
มีใครอยากลองอธิบายว่า newtonSolve[f_, guess_, steps_] := NestList[ #1 - f[#1]/f'[#1] &, guess, steps] ทำงานอย่างไรใน
Mathematica ไหมครับ เป็นการฝึกความเข้าใจเรื่อง Pure function และ Functional programming ครับ
******