一、樣條

樣條（Spline）函式是由舍恩伯格于1946年提出的，樣條是富有彈性的細木潭訓有機玻璃條，它的作用相當于“萬能”曲線板，早期船舶、汽車、飛機放樣時用鉛壓鐵壓住樣條，使其通過一系列型值點，調整壓鐵達到設計要求后繪制其曲線，稱為樣條曲線，這樣設計曲線的方法在20世紀六七十年代得到了廣泛應用，

二、幾何連續性

2.1 連續性條件

通常單一的曲線段或曲面片難以表達復雜的形狀，必須將一些曲線段拼接成組合曲線，或將一些曲面片拼接成組合曲面，才能表達復雜的形狀，為了保證在結合點處光滑過渡，需要滿足連續性條件，連續性條件有兩種：引數連續性（Parametric Continuity）與幾何連續性（Geometric Continuity），

2.2 引數連續性

• 零階引數連續性，記為C0，指相鄰兩段曲線在結合點處具有相同的坐標，
  

• 一階引數連續性，記為C1，指相鄰兩段曲線在結合點處具有相同的一階導數
  

• 二階引數連續性，記為C2，指相鄰兩段曲線在結合點處具有相同的一階導數和二階導數，
  

2.3 幾何連續性

與引數連續性不同的是，幾何連續性只要求引數成比例而非相等，

• 零階幾何連續性，記為G0，與零階引數連續性相同，即相鄰兩段曲線在結合點處有相同的坐標，
• 一階幾何連續性，記為G1，指相鄰兩段曲線在結合點處的一階導數成比例，但大小不一定相等，
• 二階幾何連續性，記為G2，指相鄰兩段曲線在結合點處的一階導數和二階導數成比例，即曲率一致，但大小不一定相等，

三、三次樣條曲線

• 三次樣條曲線是組合曲線，它在相鄰的兩型值點之間，使用三次函式進行插值，如果把n段三次曲線連接起來，使兩相鄰曲線在連接點（稱為節點）的斜率和曲率相等，就獲得n段三次函陣列合成的曲線，即三次樣條曲線，
• 三次樣條曲線是插值曲線，它通過所有型值點，
• 樣條函式的理論和運用是從三次樣條函式發展起來的，
• 在計算幾何中，應用得最早、研究得最詳細的也是三次樣條函式，因為 三次樣條函式在節點處具有C2連續性，而二階連續性是大多數工程問題所需要的，
• 樣條函式也是放樣工藝中繪制曲線用的細木條的數學模型的線性近似，符合傳統的光順要求，

3.1 三次樣條函式的定義

已知n個型值點Pi(xi, yi), i = 1, 2,…, n且a = x1 < x2 <…< xn = b，若 y=s(x)滿足下列條件：
（1）型值點Pi在函式 y=s(x) 上，
（2） s(x) 在整個區間[a, b]上二階連續可導，
（3）在每個子區間[xi, xi+1], i = 1, 2,…, n-1上，分段函式 都是引數x的三次多項式，則稱函式 是過型值點的三次樣條函式，由三次樣條函式構成的曲線稱為三次樣條曲線，

3.2 三次樣條函式的運算式

第i段的 si(x) 表示為:

式中，ai，bi，ci，di為待定系數，i=1，2 … ，n-1，

第i段曲線的首端通過Pi（xi，yi），末端通過Pi+1（xi+1，yi+1）
Pi處的二階導數為Mi，
其中各型值點的彎矩Mi的力學解釋如下：

子區間示意圖如下：

用型值點處的二階導數Mi表示的三次B樣條函式為:

3.3 求解思路

3.4 函式推導

3.5 邊界條件

3.5 推導程序

見個人資源：三次樣條插值函式求解程序.pdf

四、演算法實作

4.1 思路

（1）讀入n個型值點且滿足其x坐標遞增，
（2）根據實際情況確定三次樣條曲線的邊界條件，
（3）計算曲線的系數，將其表達為型值點二階導數的函式，
（4）用追趕法求解三彎矩方程，
（5）將求解出的引數代入三次樣條函式的運算式中，構造三次樣條曲線，
（6）回圈訪問每個節點，在每個子區間內，按照精度要求，使用直線段連接各段內的若干等分點，即可繪制出三次樣條曲線，

4.2 代碼

4.2.1 讀取型值點

首先定義資料結構p2.cpp與頭檔案(主要如下)：

CP2::CP2(double x, double y)
{
	this->x = x;
	this->y = y;
}

之后在主函式中定義讀取函式：

// CGeometricfiguretestViewmessage handlers
void CGeometricfiguretestView::ReadPoint()
{
	P[1].x = -340, P[1].y = -200;
	P[2].x = -150, P[2].y =  0;
	P[3].x =  0,   P[3].y =  -50;
	P[4].x =  100, P[4].y =  -100;
	P[5].x =  250, P[5].y =  -100;
	P[6].x =  350, P[6].y =  -50;
}

4.2.2 繪制型值點

void CGeometricfiguretestView::DrawDataPoint(CDC* pDC)//繪制型值點
{
	CBrush NewBrush, *pOldBrush;
	NewBrush.CreateSolidBrush(RGB(0, 0, 0));
	pOldBrush = pDC->SelectObject(&NewBrush);
	for( int i = 1; i < 7; i++)
		pDC->Ellipse(ROUND(P[i].x - 5), ROUND(P[i].y - 5), ROUND(P[i].x + 5), ROUND(P[i].y + 5));
	pDC->SelectObject(pOldBrush);
}

4.2.3 三次樣條曲線繪制

具體程序見注釋：

void CGeometricfiguretestView::DrawCubicSpline(CDC* pDC)//三次樣條曲線
{
	int n = 6;
	const int dim = 7;//二維陣列維數
	double b1 = 10, bn = 10;//邊界條件:"夾持端",給出起點和終點的一階導數
	double h[dim], lambda[dim], mu[dim], D[dim];
    double l[dim], m[dim], u[dim]; 
    double M[dim], K[dim];
    double a[dim], b[dim], c[dim], d[dim];
	for(int i = 1; i < n; i++)//計算hi＝xi+1-xi
    	h[i] = P[i+1].x - P[i].x;
	for(int i = 2; i < n; i++)
	{
		lambda[i] = h[i-1]/(h[i-1]+h[i]);//計算λ
		mu[i] = h[i]/(h[i-1]+h[i]);//計算μ
		D[i] = 6/(h[i-1]+h[i])*((P[i+1].y-P[i].y)/h[i]-(P[i].y-P[i-1].y)/h[i-1]);//計算D
    }
	D[1]=6*((P[2].y-P[1].y)/h[1]-b1)/h[1];//夾持端的D[1]
    D[n]=6*(bn-(P[n].y-P[n-1].y)/h[n-1])/h[n-1];//夾持端的D[n]
    mu[1]=1;//夾持端的μ[1],
	lambda[n]=1;//夾持端的λ[n]
	//追趕法求解三彎矩方程
    l[1]=2;
	u[1]=mu[1]/l[1];
    for(int i=2; i <= n; i++)
    {
		m[i]=lambda[i];
        l[i]=2-m[i]*u[i-1];
        u[i]=mu[i]/l[i];
    }
    K[1] = D[1]/l[1];//解LK=D
    for(int i = 2; i <= n;i++)
    {
		K[i]=(D[i]-m[i]*K[i-1])/l[i];
	}
	M[n] = K[n];//解UM=K
	for(int i = n-1; i >= 1;i--)
	{
		M[i]=K[i]-u[i]*M[i+1];
	}
	//計算三次樣條函式的系數
	for(int i = 1; i < n; i++)
    {
		a[i] = P[i].y;
        b[i] = (P[i+1].y-P[i].y)/h[i] - h[i]*(M[i]/3+M[i+1]/6);
        c[i] = M[i]/2;
        d[i] = (M[i+1]-M[i])/(6*h[i]);
    }
	pDC->MoveTo(ROUND(P[1].x), ROUND(P[1].y));
	double xStep = 0.5;//x的步長
	double x, y;//當前點
	for(int i = 1; i < n; i++)//回圈訪問每個節點
    {
		for(x = P[i].x; x < P[i+1].x; x += xStep)//回圈訪問每個節點
		{
			y=a[i]+b[i]*(x-P[i].x)+c[i]*(x-P[i].x)*(x-P[i].x)+d[i]*(x-P[i].x)*(x-P[i].x)*(x-P[i].x);
			pDC->LineTo(ROUND(x), ROUND(y));//繪制樣條曲線
        }
    }
}

4.2.4 主函式呼叫

void CGeometricfiguretestView::OnDraw(CDC* pDC)
{
	CTestDoc* pDoc = GetDocument();
	ASSERT_VALID(pDoc);
	if (!pDoc)
		return;
	// TODO: add draw code for native data here
	CRect rect;//定義客戶區矩形
	GetClientRect(&rect);//獲得客戶區矩形的資訊
	pDC->SetMapMode(MM_ANISOTROPIC);//自定義二維坐標系
	pDC->SetWindowExt(rect.Width(), rect.Height());//設定視窗范圍
	pDC->SetViewportExt(rect.Width(), -rect.Height());//設定視區范圍，且x軸水平向右為正，y軸垂直向上為正
	pDC->SetViewportOrg(rect.Width() / 2, rect.Height() / 2);//設定客戶區中心為二維坐標系原點
	rect.OffsetRect(-rect.Width() / 2, -rect.Height() / 2);//rect矩形與客戶區重合
	ReadPoint();
	DrawDataPoint(pDC);
	DrawCubicSpline(pDC);
}

編譯運行，可見如下：