0 目的
b樣條曲線最早是為汽車、船舶等工業產品做外觀造型設計而發明的,但是由於它有一些非常好的性質(光滑、能靈活地改變曲線形狀、表達精確),所以被用到越來越多的問題上,比如求微分方程的數值解。我們都知道,積分比微分難,對於一些複雜的問題,只會用普通的數值積分(比如尤拉積分、龍格庫塔積分法)是不夠的。一般來說,想求解稍微複雜一點的微分方程問題時(例如ivp、bvp、偏微分方程),你發現簡單的打靶法(shooting method)不好用,因為方程可能對初值非常敏感(例如敏感到需要精確到小數點後20位),所以經常採用配置法(collocation method)。要使用配置法就得選擇樣條曲線方程並需要計算曲線的導數。下面談談如何計算b樣條曲線的導數。
1 基礎知識
控制點決定了b樣條曲線的形狀,給定控制點和次數(degree),可以構造一條唯一的b樣條曲線,例如由以下6
66個控制點生成的3
33次b樣條曲線如下圖所示。如果我們改變控制點,曲線的形狀也隨之改變。具體來說,這條曲線是乙個二維的曲線(曲線上的點有 x
yx y
xy兩個座標),它也是個引數曲線,給定乙個引數值可以得到乙個二維點,把所有引數下的函式值畫出來就是圖中看到的黑色曲線。
假如你想知道在某個引數座標 u
uu 處的二維點座標應該如何計算呢?最高效的計算方法不是直接用b樣條曲線的定義式,而是採用cox-de boor公式,它是乙個遞迴形式的公式,從最底層的控制點開始經過 p
pp 步構造出最終的點, p
pp 是曲線次數。3
33次b樣條曲線(third-degree (cubic) b-spline)是最常用的b樣條曲線,這可能是因為大多數工程和物理問題都是用二階光滑微分方程描述的,選擇3
33次樣條足矣。3
33次b樣條曲線的方程總是3
33次多項式。
2 b樣條的導數
b樣條曲線不僅連續而且光滑,所以存在導數,如何求導呢?b樣條曲線的導數仍然是乙個b樣條曲線,只不過次數降低一次,控制點減少乙個。所以仍然可以像以前計算b樣條曲線一樣計算b樣條曲線的導數。而且,既然b樣條曲線的導數還是b樣條曲線,利用這個遞迴關係,更高階的導數依然是b樣條曲線。下面的例子展示了b樣條曲線的一階導數和二階導數,你可以把它們看成速度和加速度。圖中的紅色箭頭表示一階導數,黑色箭頭表示二階導數,為了方便展示,這裡只畫出了導數的方向,實際不同點處的導數幅值(箭頭長度)是不一樣的。當然還可以很容易地繼續畫出更高階的導數,這裡就不展示了。
3 解微分方程
我們來用b樣條解幾個微分方程看看效果怎麼樣。
3.1 線性常微分方程bvp
先試試簡單點的線性常微分方程,我們來解邊界值問題(bvp),即給定兩個端點的函式值。 y′′
+y′−
6y=x
x∈[0
,1]y
(0)=
0,y(
1)=1
y''+y'-6y=x \ \ \ \ \ \ \ \ \ \ x\in[0,1] \newline y(0)=0,\ \ y(1)=1
y′′+y′
−6y=
xx∈[
0,1]
y(0)
=0,y
(1)=
1這個方程存在解析解,如下[1]
^[1]y(
x)=(
43−e2
)e−3
x−(43
−1e3
)e2x
36(1e
3−e2
)−x6
−136y(x)= \frac-\left(43-\frac\right) e^}-e^2\right)}-\frac-\frac
y(x)=3
6(e3
1−e
2)(4
3−e2
)e−3
x−(4
3−e3
1)e
2x−
6x−
361
b樣條的數值解如下圖(黃色),作為對比,也畫出了解析解(黑色),相對誤差最大值為3.72797×1
0−73.72797\times10^
3.7279
7×10
−7,求解中使用的引數為collocation point數量為101,控制點數量為100。可以看到數值求解的精度還是很高的。
3.2 非線性常微分方程bvp
再試試稍微複雜點的非線性常微分方程,同樣考慮乙個邊界值問題(bvp),即被稱為bratu』s problem的 y′′
+λey
=0x∈
[0,1
]y(0
)=0,
y(1)
=0y''+ \lambda e^=0 \ \ \ \ \ \ \ \ \ \ x\in[0,1] \newline y(0)=0,\ \ y(1)=0
y′′+λe
y=0x
∈[0,
1]y(
0)=0
,y(1
)=0這個方程的解析解[2]
^[2]如下y (x
)=−2
ln[cosh((
x−12
)θ2)
cosh(θ
4)]y(x)=-2\text\left[\frac((x-\frac)\frac )}(\frac)}\right]
y(x)=−
2ln[cosh(4
θ)cosh((
x−21
)2θ
)]
引數λ
=3.513830719
\lambda=3.513830719
λ=3.51
3830
719、θ
=4.79871456
\theta=4.79871456
θ=4.79
8714
56時的b樣條的數值解如下圖(黃色),作為對比,也畫出了解析解(黑色),相對誤差最大值為2.77846×1
0−32.77846\times10^
2.7784
6×10
−3,求解中使用的引數為collocation point數量為101,控制點數量為100。
B 樣條曲線的導數
回目錄 定義如下 如果原始的clamped節點向量是u 0 p 1 u p 1,u m p 1,u m p 1 那麼移動第乙個和最後乙個節點使得第乙個和最後乙個節點重複度變成,p 而不是p 1,我們有乙個m 1 個節點u 0 p u p 1,u m p 1,u m p 的新節點序列。那麼,可證明在原...
樣條曲線 B樣條曲線求解及C 實現
功能 根據引數u值和k 大小為階數值 與節點向量,計算第i個k次b樣條基數 輸入引數 u 引數值 k 大小值為階數 i 第i個k次b樣條的支撐區間左端節點的下標 anode為節點向量。輸出引數 返回函式值。double getbasefunval double u,int i,int k,vecto...
b樣條和三次樣條 樣條曲線
最近在學習軌跡規劃中的軌跡生成,涉及到樣條曲線方面的知識,總結一下。曲線的平滑性和相應的平滑性的評判準則相關,在 1 中,作者採用曲率的平方和曲率導數的平方作為評判準則 其中 是路徑點的方向角。最小化這兩個準則的軌跡分別是圓弧和三階螺旋線,並對在對稱和不對稱情況下如何生成路徑進行了分析,事實表明三階...