計算物理基礎(chǔ)-解方程(組)
? ? ? ?在解決物理問題時,我們通常引入一定數(shù)目的變量,根據(jù)已知條件建立方程、方程組,然后通過數(shù)學(xué)方法來求解。除了有求根公式的方程、線性方程組等,在多數(shù)情況下我們需要借助數(shù)值方法通過自洽迭代來處理,這里給出常用的庫函數(shù)和迭代程序范例供參考。
(1)庫函數(shù)roots和fsolve
? ? ? ?對于多項式求根,先初始化系數(shù)向量c,然后用roots函數(shù)即可。這里求解域是復(fù)數(shù),對于x^3-x-1=0有3個根,其中一個是實數(shù)根,虛部為零。

? ? ?也可以用fsolve函數(shù),對應(yīng)的語法是fun = @(x) sqrt(x.^3-x-1); fsolve(fun,1),返回結(jié)果為?1.3247,只有實數(shù)根。注意fsolve(fun,1)這里的1是初值,如果用fsolve(fun,0.5)將返回結(jié)果?-0.5774,并不是方程的根,因為迭代過程和初始值是相關(guān)的。

? ? 這里p_1,p_2,x_1,x_2是4個未知數(shù),對應(yīng)x向量中4個分量,[rand,rand,rand,rand]是四個變量的初始值取了(0,1)之間的隨機(jī)數(shù),計算結(jié)果為?1.0000? ?1.0000? -0.5773? ?0.5773,是之前Gauss求積法里面待定系數(shù)的數(shù)值。
(2)牛頓法、割線法
? ? ? ?假定在x(k+1)時得到方程的解,泰勒展開得到x(k+1)和x(k)的遞推關(guān)系,注意這里需要用到f(x)的導(dǎo)數(shù)。對于非線性方程組,此時導(dǎo)數(shù)用雅可比矩陣代替,除以導(dǎo)數(shù),變成求雅可比矩陣的逆。在牛頓法中,用差分形式代替導(dǎo)數(shù),對應(yīng)的是割線法,需要用兩個初始值。

(3)復(fù)數(shù)域的求根
? ? ? ?將函數(shù)拓展到復(fù)數(shù)域,限制實部和虛部范圍都在(-2,2)之間。從任意的一點出發(fā),通過牛頓法迭代,將演化到3個復(fù)數(shù)根之一,分別用不同的顏色表示。值得注意的是,初始值分布無法用簡單的區(qū)域分割,在某些區(qū)域,兩個相差很小的初值,經(jīng)過迭代,可能收斂到不同的復(fù)數(shù)根,在平面上的分布有類似分形的結(jié)構(gòu)。

附錄:
A. 牛頓法解方程
clear
x(1)=-10;x(2)=x(1)-(x(1)^3-x(1)-1)/(3*x(1)*x(1)-1);
i=1;
while abs(x(i+1)-x(i))>0.001?
i=i+1;?
x(i+1)=x(i)-(x(i)^3-x(i)-1)/(3*x(i)*x(i)-1);
end
plot(x,'o-')
B. 割線法解方程
clear?
x(1)=1;x(2)=2;i=2;
while abs(x(i)-x(i-1))>0.01
f(i)=x(i)^3-x(i)-1;
f(i-1)=x(i-1)^3-x(i-1)-1;
x(i+1)=x(i)-(f(i)*(x(i)-x(i-1)))/(f(i)-f(i-1));?
i=i+1;?
end