一.非線性方程組
在前面我們討論的大多數是關於未知數的線性方程組，本節主要討論非線性方程組，簡單的理解為不是關於未知數（不一定是x）的一次方程，比如我們熟悉的x^2-4x+3=0這種方程組但是這種方程比較簡單，我們可以利用一些公式直接解出其解，但是對於x^3-sin(x)=0這種方程，我們是沒有通用公式的解出它的解，這時我們就要採取一些其他方法來求解。本節我們主要採取的是採用迭代循環方法進行求解，比較典型方法有二分法和牛頓法。(其實這兩種方法的基礎就是通過畫圖來大概確定解範圍，從而在相應區間進行求解）
一.畫圖求解範圍
import numpy as npimport matplotlib.pyplot as pltx=np.linspace(-2,2,1000)
f1=x**3-x**2-xf2=x-3 * np.sin(x)f3=np.exp(x)-2f4=np.sin(x)-x**2
fig,axes=plt.subplots(1,4,figsize=(12, 3),sharey=True)
for n,f in enumerate([f1,f2,f3,f4]): axes[n].plot(x,f,lw=1.5) axes[n].axhline(0,ls=':',color='k') axes[n].set_ylim(-5,5) axes[n].set_xticks([-2,-1,0,1,2]) axes[n].set_xlabel('x',fontsize=18)
axes[0].set_ylabel('f(x)',fontsize=18)
titles = [r'$f(x)=x^3-x^2-x$',r'$f(x)=x-3*sin(x)$', r'$f(x)=\exp(x)-2$',r'$f(x)=sin(x)-x^2$']for n,title in enumerate(titles): axes[n].set_title(title) fig.tight_layout()plt.show()運行結果
二.二分法求解
import numpy as npimport matplotlib.pyplot as pltplt.rcParams["font.family"] = "SimHei"plt.rcParams["font.size"] = "12"f=lambda x: np.exp(x) - 2tol=0.1a,b=-2,2x=np.linspace(-2.1, 2.1, 1000)
fig,ax=plt.subplots(1,1,figsize=(12, 4))
ax.plot(x,f(x),lw=1.5)
ax.axhline(0,ls=':',color='k')ax.set_xticks([-2,-1,0,1,2])ax.set_xlabel('x',fontsize=18)ax.set_ylabel('f(x)',fontsize=18)
fa,fb=f(a),f(b)
ax.plot(a,fa,'ko')ax.plot(b,fb,'ko')ax.text(a,fa+0.5,"a",ha='center',fontsize=18)ax.text(b,fb+0.5,"b",ha='center',fontsize=18)
n=1while b-a>tol: m=a+(b-a)/2 fm=f(m) ax.plot(m,fm,'ko') ax.text(m,fm-0.5,r"$m_%d$"%n,ha='center') n+=1 if np.sign(fa)==np.sign(fm): a,fa=m,fm else: b,fb=m,fm
ax.plot(m,fm,'r.',markersize=10)ax.annotate("方程根的大約範圍 %.3f"% m, fontsize=14,family="SimHei", xy=(a,fm),xycoords='data', xytext=(-150,+50),textcoords='offset points', arrowprops=dict(arrowstyle="->",connectionstyle="arc3,rad=-.5"))
ax.set_title("二分法求解方程組")plt.show()print('總的迭代次數為%d'%n)作者對原始碼進行了改進
運行結果
總的迭代次數為7（由於我們之前已經畫圖大概確定了它的解範圍，因此迭代效率總體來說還是很高的）
三.牛頓法
import numpy as npimport matplotlib.pyplot as pltimport sympyplt.rcParams["font.family"] = "SimHei"plt.rcParams["font.size"] = "12"
tol=0.01xk=2
s_x=sympy.symbols("x")s_f=sympy.exp(s_x)-2
f=lambda x:sympy.lambdify(s_x,s_f)(x)fp=lambda x:sympy.lambdify(s_x,sympy.diff(s_f,s_x))(x)
x=np.linspace(-1,2.1,1000)
fig,ax=plt.subplots(1,1,figsize=(12,4))ax.plot(x,f(x))
ax.axhline(0,ls=':',color='k')
n=0while f(xk)>tol: xk_new=xk-f(xk)/fp(xk) ax.plot([xk,xk],[0,f(xk)],color='k',ls=':') ax.plot(xk,f(xk),'ko') ax.text(xk,-.5,r'$x_%d$'% n,ha='center') ax.plot([xk,xk_new],[f(xk),0],'k-') xk=xk_new n+=1
ax.plot(xk,f(xk),'r.',markersize=15)ax.annotate("方程根的大約範圍為%.3f"%xk, fontsize=14,family="SimHei", xy=(xk,f(xk)),xycoords='data', xytext=(-150,+50),textcoords='offset points', arrowprops=dict(arrowstyle="->",connectionstyle="arc3,rad=-.5"))
ax.set_title("牛頓法求解")ax.set_xticks([-1,0,1,2])fig.tight_layout()plt.show()print('總的迭代次數為%d'%n)作者對原始碼進行了改進
運行結果
總的迭代次數為4（由於我們之前已經畫圖大概確定了它的解範圍，因此迭代效率總體來說還是很高的）