iT邦幫忙

第 11 屆 iT 邦幫忙鐵人賽

DAY 7
0
AI & Data

連接數學與現實世界的橋樑 -- 數學建模系列 第 7

Day 7 : 當方程組無法用手算怎麼辦 -- 數值分析

在此之前,我們所舉的例子都可以透過人工手算的方式求解,但是當方程式很複雜該怎麼辦,如
https://chart.googleapis.com/chart?cht=tx&chl=%5Cleft%20%5C%7B%2015.5x%5E%7B-0.5%7D%20-%208%20%2B%2013y%5E%7B-0.2%7D%20-%200.064yx%5E%7B-1.08%7D%20%3D%200%20%5C%5C%209y%5E%7B-0.4%7D%20-%205%20%2B%200.8x%5E%7B-0.08%7D%20-%200.26xy%5E%7B-1.2%7D%20%3D%200
這時可以透過數值分析來幫助我們求解,這裡我們無法去詳述所有的數值分析理論與方法,所以在此就大概介紹一下兩種方法,牛頓法與隨機搜尋法。

目標函數:
https://chart.googleapis.com/chart?cht=tx&chl=200ce%5E%7Bcx%7D(0.65%20-%200.01x)%20-%202e%5E%7Bcx%7D%20-%200.45,令https://chart.googleapis.com/chart?cht=tx&chl=c%20%3D%200.025,求x,這樣複雜的函數很難用手算,求得一個解析解。此時,可以透過牛頓法來求得https://chart.googleapis.com/chart?cht=tx&chl=x的值。

牛頓法

簡單來說就是,在解的附近,利用線性近似來逼近https://chart.googleapis.com/chart?cht=tx&chl=F(x),亦即https://chart.googleapis.com/chart?cht=tx&chl=F(x)%20%5Capprox%20F(x_0)%20%2B%20F'(x_0)(x%20-%20x_0)。此時,https://chart.googleapis.com/chart?cht=tx&chl=x_0就是函數的近似解。
為了得到更好的近似,就會透過迭代的方式,找出更靠近解的近似值。作法如下,
https://chart.googleapis.com/chart?cht=tx&chl=F(x_0)%20%2B%20F'(x_0)(x%20-%20x_0)%20%3D%200
解出https://chart.googleapis.com/chart?cht=tx&chl=x%20%3D%20x_1%20%20%5CRightarrow%20%20x_1%20%3D%20x_0%20-%20%5Cfrac%20%7BF(x_0)%7D%7BF'(x_0)%7D

變量:
x(n) = n次迭代後的根的近似解
N = 迭代次數
Pseudo code:
For n = 1~ N:
https://chart.googleapis.com/chart?cht=tx&chl=%5Cqquad%20x(n)%20%5Cleftarrow%20x(n-1)%20-%20F(x(n-1))%20%2F%20F'(x(n-1))
end
Python

# 目標函數
def func(c, x):
    return 200*c*math.exp(c*x)*(0.65 - 0.01*x) - 2*math.exp(c*x)-0.45

# 目標函數的一次導數
def der_func(c, x, h=0.01):
    return (func(c, x+h) - func(c, x)) / h

N = 100               # 迭代次數
initial_param = 0.0   # 初始值
final_param = 0.0     # 終值
diff = 10e5           # 差異
c=0.025               # 參數
cnt = 0               # 達標次數
error = []            # 迭代過程差異的紀錄

while (abs(diff) >= 10e-5) and (cnt <= N):
    f = func(c, initial_param)
    diff_f = der_func(c, initial_param)
    final_param = initial_param - (f/diff_f)
    diff = final_param - initial_param
    error.append(abs(diff))
    initial_param = final_param
    cnt += 1
    
print(final_param, cnt)

隨機搜尋法

這方法比起牛頓法更加簡單,就是從解空間中隨機抽取N個點,去看哪個點離極值最近。此法的準確度跟抽取的點數個數有關,抽樣次數越多,越可能找到極值,但相對的也越花時間。
Python

def func(c, x):
    return 200*c*math.exp(c*x)*(0.65 - 0.01*x) - 2*math.exp(c*x)-0.45

min_y = 10e5
N = 1000
threshold = 10e-5
c=0.025
final_x = 0

for i in range(N):
    approx_x = random.uniform(10, 30)
    approx_y = func(c, approx_x)
    
    if abs(approx_y) < min_y:
        min_y = abs(approx_y)
        final_x = approx_x
    
print(final_x, min_y)

這裡只介紹了兩種單變量的數值方法,關於多變量及更多的數值求解法可以參考數值方法的教科書。


上一篇
Day 6 : 有限制條件之最佳化 - KKT定理
下一篇
Day 8 : 全域最佳化 - 凸函數
系列文
連接數學與現實世界的橋樑 -- 數學建模30

尚未有邦友留言

立即登入留言