微分方程的数值解法与程序实现 pdf_数值计算方法·第三部分
第三部分內容主要是非線性的各種玩意兒,比如一個非線性的方程怎么解、一個非線性的函數曲線怎么擬合、一個積分怎么算、一個微分方程怎么求數值解:
首先一定要判斷求根區間上是否有重根/好幾個單根,或者某個根是不是重根!觀察原來的方程
:重根可以這樣定義:針對單根區間(自變量的取值范圍里,方程只有一個根)的情況,非線性方程的解法可以參考線性方程組:既然
能變換出 這種形式, 也可以化成 ,而這個方程迭代下去,最終如果取得 ,那么 就是原方程的一個根了理想的情況下這種迭代法并不一定收斂,要求是:
也就是要求在單根區間上, 的值在 里取, 的值在 里取不滿足(1)的情況不滿足(2)的情況附加了一個誤差的結果,編程算的時候可以判斷什么時候停止迭代。作為了解
需要注意,這樣即使證明收斂,也是在單根區間上局部收斂,別處收不收斂還是不清楚的
不如破罐子破摔吧!直接判斷是否在根
附近是不是收斂,有個巴掌大地方收斂就夠了,順便用收斂的階來判斷下收斂速度: 常數,則稱p階收斂還有種簡單的判斷方法: ,線性收斂(一階收斂), 越接近0,收斂越快 ,二階收斂
……(必須p階連續可導,前p-1階導數都是0),p階收斂
求到哪一階的導數,代入 后不是0了,說明哪一階收斂
最常用的是Newton法,這個至少二階收斂:
(單根區間上的牛頓法/切線法)由于Newton法有時候求導數有點麻煩,直接用差商
代替導數往里代,或者還用直接用初始點 構造公式的,都作為了解即可:Newton法局部收斂,初始點選不好容易發散,也可以在被減的那一項前面加個系數
,保證每一步函數值 逐步接近0,也作為了解:每一次迭代都要重新算λ對多根區間的處理,一是這個區間上所有根互不相同,僅有單根,二是有重根:
- 對只有單根的情況,可以分出跟單根數目一樣多的單根區間:
在單根區間上求根,用二分法就可以:
最后剩下的區間足夠小的時候,取中點作根的近似值就行- 對重根,此時用Newton法收斂階數下降了,也就是收斂更慢了,在被減項前加一個根的重數 ,還能按二階收斂的速度解出根:
一個非線性的函數曲線怎么擬合呢?有兩種想法比較簡單:
1.把這些點連起來,畫出一條曲線——這叫插值
2.不連線了,直接在散點中間隨便畫條線,只要所有的點都落在線不遠的地方就行——這叫逼近
下面假設已知函數的 個點,從0開始編號,這樣得到的插值和逼近函數都是 次的對插值,還分為在對應位置取對應的函數值,但和原來的函數導數可以不一樣(如Lagrange插值和Newton插值)、導數也跟原函數一樣(Hermite插值)、導數盡量連續(這樣擬合出的曲線光溜溜的沒有棱角,如樣條插值)
- Lagrange插值和Newton插值(剩下的插值方法可以看成其變種)
在掌握這兩個插值方法之前,首先要確定什么時候可以做這種插值,什么時候不可以:
我們要求插值所得的函數(插值函數)在節點上和被插值函數相等,即
也就是找一組函數(插值基函數)
的線性組合 ,讓它在節點上等于被插值函數對兩種插值來說,
,所以實際上就是求 ,節點數應該跟 的數目一樣多,就相當于解這樣一個方程組:(Haar條件)這樣一個方程組有解的條件是:在確認插值可以做的情況下就可以弄出插值的公式了,兩個插值實際上相同,只不過用于計算的公式不同罷了:
Lagrange的插值函數是式中
(連乘, 分子去掉了 ,分母去掉了 )Newton的插值函數是
式中
(k階差商的定義,中間 任意調換順序都行)
對Newton插值,可以記住這樣一張均差表:
兩個插值函數對原函數來說都不是完全精確的(開局幾個點,誰?知道到底是哪條線呢?),誤差也可以估計出來:
作了解就行,知道n次的插值多項式能夠精確表示≤n次的函數,也就是n階精度也就是說,一個 次的插值多項式如果用于插值 次及以下的函數,和原函數是完全相等的用于應對這種問題很方便,前一個空是0,后一個空是x^3+1兩個插值(實際上得出同一個插值函數)都沒有很好的收斂性,會出現龍格(Runge)現象:
Runge現象:插值節點一多,次數變高,插值函數兩頭開始抽搐了起來~所以要避免插值節點變多,可以采取的好方法是中間分成好些段,在這些段上只用兩三個點來做插值,段和段之間井水不犯河水。作為了解即可
Hermite插值是上述兩種插值的變種,是在已知被插值函數的一些點的導數(二階導、三階導的情況下),把這些信息考慮進插值函數而作的一種調整,這樣圖像也在一定程度上變光滑了,可以用以下的均差表來算:
同樣作為了解即可三次樣條插值利用次數不超過三次的樣條函數作為插值函數,無需知道被插值函數的導數,樣條函數是一種分段函數,它的前后兩段2階連續可導我們要學會判斷什么是三次樣條函數:
能使前后兩段2階連續可導,就要在兩段的交點
處有 ,我們把后一段 寫成: ,也就有 我們用每一段的表達式減前一段的表達式,如果減出來都是上面那個形式,就可以證明這是一個三次樣條函數剩下的作為了解即可
- (標準)正交多項式和最小二乘法逼近
我們定義內積
,其中權函數是 ,它可積、非負且零點只有有限個如果上面那個內積為0,則稱 和 在 上正交(正交多項式一定要講在哪個區間上,以什么作為權函數)如果 則稱 是標準的
設有一組線性無關的函數
,通過以下方法得到正交函數 :(把最后一列替換成函數,如果 ,看這個行列式的右下角就知道函數的次數)
最常見的是冪函數
生成的正交函數系,它們具有這樣的性質: 是 次多項式,與低于 次的所有多項式(不僅是正交多項式!)正交,并且在 里有 個互異零點遞推公式性質作為了解 上權函數是 的正交多項式稱為Chebyshev多項式,具有以下性質:零點公式可以記住,計算挺快的經過正交化還可以得到標準正交函數系,作為了解:
最常見的最小二乘法是利用冪函數
作為基,求解 中 的方法,要求所有離散點的 的和最小(“最小二乘”),因為是離散點,定義 為離散型內積(連續型內積是積分,離散型內積是求和) 稱為法方程組,實際解這個方程組即可其中
這種形式挺重要的,對應上式化簡成:對非線性的函數,只要化成上述的線性形式,把自變量和因變量轉換過來,就可以利用上式化簡
知道怎么用一個好算的函數代替一個非線性的函數,那我們用這個好算的函數來算后者的微分和積分,不也很容易嗎?
課程主要講的是數值積分的內容,微分和積分的思想差不多,就是用有限個點的函數值湊成個公式,代替原來的整段函數去算
確定
和 大概有兩種方法:靠插值(一般取等距節點)、靠逼近(主要用的是上一章的正交多項式)在介紹具體方法之前,先要理解代數精度:
到底在什么情況上可以劃等號多項式函數并不是世界的全部,但是它挺簡單,適合去衡量一種方法- 插值-等距節點
直接把積分
里面那個函數換成插值函數就行,然后就能導出下面的結果(Newton-Cotes公式):了解就行直接確定
中每個 的值,就有了下面三個公式:也不用記,直接用復化求積公式n=1的情況就行三者的代數精度是1,3,5上述等距節點求積公式用的是插值函數,復化求積公式用的則是分段低次插值函數,在每個子區間上套用上述結果,然后共有
個區間,對應就能導出:這三個公式的特點是:分母 是多少,就有多少個 相加;在給定有限點的情況下,在 的公式里算成 的點,在 里算 ,在 里算
- 逼近-Gauss點
所以求積節點只要求出正交多項式 ,然后求出它的根就行
而其求積系數則令 ,代入 計算求積系數性質了解一下,知道總為正數即可
Gauss-Chebyshev求積公式(區間
,權函數 )常用,記住以下結論:n=點數-1常用的公式是
(兩點,代數精度為3)和 (三點,代數精度為5)的:不在
區間上的積分想借助這些結論去算,可以作變量替換:最后就是用數值方法解微分方程了。
研究一階常微分方程
,通常我們知道這個方程在一個點上的解 ,設置一個步長 ,這樣一步一步往下探,解出 的過程。設在對應的 的位置上,精確的解應該是常用的迭代方法可以寫成以下的通式:
(這個方程就是用來求 的,并且已知條件從 或 開始寫起,所以有后面的兩個限制條件)
從這個方程可以看出:
1.迭代是“幾步法”: 是幾(需要前面的幾個值)
2.顯式(不需要通過解方程求得 )和隱式: 即為顯性的,反之為隱形的
在了解具體方法之前,同樣有一個誤差的概念:
假設前
步沒有誤差, 稱為局部截斷誤差前
步誤差累積, 稱為整體截斷誤差假設 (原函數是一個次數低于 的多項式)局部截斷誤差可以寫成多項式 的形式
這兩項當中的 稱為局部截斷誤差主項, 是局部截斷誤差主項系數
整體截斷誤差比局部截斷誤差低一階,稱對應的迭代方法為 階的
用Taylor展開法展開上面的通式,可以得到:
找到下標最大的那個等于0的 ,它的下標就是階數
除了誤差之外,收斂性(步長趨近于0時是否求得精確值)、絕對穩定性(計算時的舍入誤差是否影響特別大)也很重要,有以下結論:
(1) 的所有根都在單位圓內或圓上,且圓上的根是單根;(2)階數至少為1滿足以上兩個條件,迭代方法是收斂的(常用的方法都是收斂的)用模型問題 的模型問題判斷是否絕對穩定,并令實部為負的復數 ,能令誤差不擴大的 取值范圍就是絕對穩定區域,它與實軸的交集就是絕對穩定區間
保證 的所有根都在單位圓內,利用下述方法得出 的取值范圍,就能解出絕對穩定區域和絕對穩定區間:把方程化成上述形式,找出b,c,注意是減號
常用的幾種單步法可以記住它的階數和絕對穩定區間:
顯式Euler法: ,階數為1,絕對穩定區間隱式Euler法: (需要解方程),階數為1,絕對穩定區間
梯形法: (需要解方程),階數為2,絕對穩定區間 改進的Euler法,實際上把f_{n+1}拆開了
總結
以上是生活随笔為你收集整理的微分方程的数值解法与程序实现 pdf_数值计算方法·第三部分的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 前端学习(2592):当前用户显示
- 下一篇: 工作177:表单重置项目处理