情報学部 菅沼ホーム 目次 索引

補間法（スプライン）

    1. A. C++
    2. B. Java
    3. C. JavaScript
    4. D. PHP
    5. E. Ruby
    6. F. Python
    7. G. C#
    8. H. VB

　　プログラムは，3 次スプライン関数によってスプライン補間するためのものです．なお，このプログラムは，クラスを使用して記述しています（ JavaScript の場合はオブジェクト）．

  1. C++ 

    /********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
#include <stdio.h>

double spline(int, double *, double *, double,
              double *, double *, double *, double *, double *, double *);

int main()
{
	double *h, *b, *d, *g, *u, *r, *x, *y, sp, x1, y1;
	int n;
					// 設定
	n = 5;
	h = new double [n];
	b = new double [n];
	d = new double [n];
	g = new double [n];
	u = new double [n];
	r = new double [n+1];
	x = new double [n+1];
	y = new double [n+1];

	x[0] = 0.0;
	x[1] = 1.0;
	x[2] = 2.0;
	x[3] = 3.0;
	x[4] = 4.0;
	x[5] = 5.0;
	y[0] = 0.0;
	y[1] = 1.0;
	y[2] = 4.0;
	y[3] = 9.0;
	y[4] = 16.0;
	y[5] = 25.0;
					// 実行と出力
	x1 = 0.0;
	sp = 0.1;
	while (x1 <= 5.01) {
		y1 = spline(n, x, y, x1, h, b, d, g, u, r);
		printf("%f %f\n", x1, y1);
		x1 += sp;
	}

	delete [] h;
	delete [] b;
	delete [] d;
	delete [] g;
	delete [] u;
	delete [] r;
	delete [] x;
	delete [] y;

	return 0;
}

/**********************************/
/* ３次スプライン関数による補間   */
/*      n : 区間の数              */
/*      x, y : 点の座標           */
/*      x1 : 補間値を求める値     */
/*      h, b, d, g, u, r : 作業域 */
/*      return : 補間値           */
/*      coded by Y.Suganuma       */
/**********************************/
double spline(int n, double *x, double *y, double x1,
              double *h, double *b, double *d, double *g, double *u, double *r)
{
	int i = -1, i1, k;
	double y1, qi, si, xx;
					// 区間の決定
	for (i1 = 1; i1 < n && i < 0; i1++) {
		if (x1 < x[i1])
			i = i1 - 1;
	}
	if (i < 0)
		i = n - 1;
					// ステップ１
	for (i1 = 0; i1 < n; i1++)
		h[i1] = x[i1+1] - x[i1];
	for (i1 = 1; i1 < n; i1++) {
		b[i1] = 2.0 * (h[i1] + h[i1-1]);
		d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1]);
	}
					// ステップ２
	g[1] = h[1] / b[1];
	for (i1 = 2; i1 < n-1; i1++)
		g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
	u[1] = d[1] / b[1];
	for (i1 = 2; i1 < n; i1++)
		u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
					// ステップ３
	k      = (i > 1) ? i : 1;
	r[0]   = 0.0;
	r[n]   = 0.0;
	r[n-1] = u[n-1];
	for (i1 = n-2; i1 >= k; i1--)
		r[i1] = u[i1] - g[i1] * r[i1+1];
					// ステップ４
	xx = x1 - x[i];
	qi = (y[i+1] - y[i]) / h[i] - h[i] * (r[i+1] + 2.0 * r[i]) / 3.0;
	si = (r[i+1] - r[i]) / (3.0 * h[i]);
	y1 = y[i] + xx * (qi + xx * (r[i]  + si * xx));

	return y1;
}
    クラスを使用した場合
    /********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
#include <stdio.h>

/****************/
/* クラスの定義 */
/****************/
class Spline
{
		double *r, *q, *s, *x, *y;
		int n;

	public:

		/**************************/
		/* コンストラクタ         */
		/*      n1 : 区間の数     */
		/*      x1, y1 : 点の座標 */
		/**************************/
		Spline(int n1, double *x1, double *y1)
		{
			int i1;
						// 領域の確保
			n = n1;

			double *h = new double [n];
			double *b = new double [n];
			double *d = new double [n];
			double *g = new double [n];
			double *u = new double [n];
			q = new double [n];
			s = new double [n];
			r = new double [n+1];
			x = new double [n+1];
			y = new double [n+1];

			for (i1 = 0; i1 <= n; i1++) {
				x[i1] = x1[i1];
				y[i1] = y1[i1];
			}
						// ステップ１
			for (i1 = 0; i1 < n; i1++)
				h[i1] = x[i1+1] - x[i1];
			for (i1 = 1; i1 < n; i1++) {
				b[i1] = 2.0 * (h[i1] + h[i1-1]);
				d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1]);
			}
						// ステップ２
			g[1] = h[1] / b[1];
			for (i1 = 2; i1 < n-1; i1++)
				g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
			u[1] = d[1] / b[1];
			for (i1 = 2; i1 < n; i1++)
				u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
						// ステップ３
			r[0]   = 0.0;
			r[n]   = 0.0;
			r[n-1] = u[n-1];
			for (i1 = n-2; i1 >= 1; i1--)
				r[i1] = u[i1] - g[i1] * r[i1+1];
						// ステップ４
			for (i1 = 0; i1 < n; i1++) {
				q[i1] = (y[i1+1] - y[i1]) / h[i1] - h[i1] * (r[i1+1] + 2.0 * r[i1]) / 3.0;
				s[i1] = (r[i1+1] - r[i1]) / (3.0 * h[i1]);
			}

			delete [] h;
			delete [] b;
			delete [] d;
			delete [] g;
			delete [] u;
		}

		/****************/
		/* デストラクタ */
		/****************/
		~Spline()
		{
			delete [] q;
			delete [] s;
			delete [] r;
			delete [] x;
			delete [] y;
		}

		/******************************/
		/* 補間値の計算               */
		/*      x1 : 補間値を求める値 */
		/*      return : 補間値       */
		/******************************/
		double value(double x1)
		{
			int i = -1, i1;
			double y1, xx;
						// 区間の決定
			for (i1 = 1; i1 < n && i < 0; i1++) {
				if (x1 < x[i1])
					i = i1 - 1;
			}
			if (i < 0)
				i = n - 1;
						// 計算
			xx = x1 - x[i];
			y1 = y[i] + xx * (q[i] + xx * (r[i]  + s[i] * xx));

			return y1;
		}
};

/*************/
/* main 関数 */
/*************/
int main()
{
	double *x, *y, sp, x1, y1;
	int n;
					// 設定
	n = 5;
	x = new double [n+1];
	y = new double [n+1];

	x[0] = 0.0;
	x[1] = 1.0;
	x[2] = 2.0;
	x[3] = 3.0;
	x[4] = 4.0;
	x[5] = 5.0;
	y[0] = 0.0;
	y[1] = 1.0;
	y[2] = 4.0;
	y[3] = 9.0;
	y[4] = 16.0;
	y[5] = 25.0;

	Spline spline(n, x, y);
					// 実行と出力
	x1 = 0.0;
	sp = 0.1;
	while (x1 <= 5.01) {
		y1 = spline.value(x1);
		printf("%f %f\n", x1, y1);
		x1 += sp;
	}

	delete [] x;
	delete [] y;

	return 0;
}

  2. Java 

    /********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
import java.io.*;

public class Test {
	public static void main(String args[]) throws IOException
	{
		double h[], b[], d[], g[], u[], r[], sp, x1, y1;
		int n;
					// 設定
		n = 5;
		h = new double [n];
		b = new double [n];
		d = new double [n];
		g = new double [n];
		u = new double [n];
		r = new double [n+1];
		double x[] = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0};
		double y[] = {0.0, 1.0, 4.0, 9.0, 16.0, 25.0};
					// 実行と出力
		x1 = 0.0;
		sp = 0.1;
		while (x1 <= 5.01) {
			y1 = spline(n, x, y, x1, h, b, d, g, u, r);
			System.out.println(x1 + " " + y1);
			x1 += sp;
		}
	}

	/**********************************/
	/* ３次スプライン関数による補間   */
	/*      n : 区間の数              */
	/*      x, y : 点の座標           */
	/*      x1 : 補間値を求める値     */
	/*      h, b, d, g, u, r : 作業域 */
	/*      return : 補間値           */
	/*      coded by Y.Suganuma       */
	/**********************************/
	static double spline(int n, double x[], double y[], double x1,
           double h[], double b[], double d[], double g[], double u[], double r[])
	{
		int i = -1, i1, k;
		double y1, qi, si, xx;
					// 区間の決定
		for (i1 = 1; i1 < n && i < 0; i1++) {
			if (x1 < x[i1])
				i = i1 - 1;
		}
		if (i < 0)
			i = n - 1;
					// ステップ１
		for (i1 = 0; i1 < n; i1++)
			h[i1] = x[i1+1] - x[i1];
		for (i1 = 1; i1 < n; i1++) {
			b[i1] = 2.0 * (h[i1] + h[i1-1]);
			d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1]);
		}
					// ステップ２
		g[1] = h[1] / b[1];
		for (i1 = 2; i1 < n-1; i1++)
			g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
		u[1] = d[1] / b[1];
		for (i1 = 2; i1 < n; i1++)
			u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
					// ステップ３
		k      = (i > 1) ? i : 1;
		r[0]   = 0.0;
		r[n]   = 0.0;
		r[n-1] = u[n-1];
		for (i1 = n-2; i1 >= k; i1--)
			r[i1] = u[i1] - g[i1] * r[i1+1];
					// ステップ４
		xx = x1 - x[i];
		qi = (y[i+1] - y[i]) / h[i] - h[i] * (r[i+1] + 2.0 * r[i]) / 3.0;
		si = (r[i+1] - r[i]) / (3.0 * h[i]);
		y1 = y[i] + xx * (qi + xx * (r[i]  + si * xx));

		return y1;
	}
}
    クラスを使用した場合
    /********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
import java.io.*;

public class Test_c {
	public static void main(String args[]) throws IOException
	{
		double sp, x1, y1;
		int n;
					// 設定
		n = 5;
		double x[] = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0};
		double y[] = {0.0, 1.0, 4.0, 9.0, 16.0, 25.0};

		Spline spline = new Spline(n, x, y);
					// 実行と出力
		x1 = 0.0;
		sp = 0.1;
		while (x1 <= 5.01) {
			y1 = spline.value(x1);
			System.out.println(x1 + " " + y1);
			x1 += sp;
		}
	}
}

/********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
class Spline
{
	private double r[], q[], s[], x[], y[];
	private int n;

	/**************************/
	/* コンストラクタ         */
	/*      n1 : 区間の数     */
	/*      x1, y1 : 点の座標 */
	/**************************/
	Spline(int n1, double x1[], double y1[])
	{
		int i1;
						// 領域の確保
		n = n1;

		double h[] = new double [n];
		double b[] = new double [n];
		double d[] = new double [n];
		double g[] = new double [n];
		double u[] = new double [n];
		q = new double [n];
		s = new double [n];
		r = new double [n+1];
		x = new double [n+1];
		y = new double [n+1];

		for (i1 = 0; i1 <= n; i1++) {
			x[i1] = x1[i1];
			y[i1] = y1[i1];
		}
						// ステップ１
		for (i1 = 0; i1 < n; i1++)
			h[i1] = x[i1+1] - x[i1];
		for (i1 = 1; i1 < n; i1++) {
			b[i1] = 2.0 * (h[i1] + h[i1-1]);
			d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1]);
		}
						// ステップ２
		g[1] = h[1] / b[1];
		for (i1 = 2; i1 < n-1; i1++)
			g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
		u[1] = d[1] / b[1];
		for (i1 = 2; i1 < n; i1++)
			u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
						// ステップ３
		r[0]   = 0.0;
		r[n]   = 0.0;
		r[n-1] = u[n-1];
		for (i1 = n-2; i1 >= 1; i1--)
			r[i1] = u[i1] - g[i1] * r[i1+1];
						// ステップ４
		for (i1 = 0; i1 < n; i1++) {
			q[i1] = (y[i1+1] - y[i1]) / h[i1] - h[i1] * (r[i1+1] + 2.0 * r[i1]) / 3.0;
			s[i1] = (r[i1+1] - r[i1]) / (3.0 * h[i1]);
		}
	}

	/******************************/
	/* 補間値の計算               */
	/*      x1 : 補間値を求める値 */
	/*      return : 補間値       */
	/******************************/
	double value(double x1)
	{
		int i = -1, i1;
		double y1, xx;
						// 区間の決定
		for (i1 = 1; i1 < n && i < 0; i1++) {
			if (x1 < x[i1])
				i = i1 - 1;
		}
		if (i < 0)
			i = n - 1;
						// 計算
		xx = x1 - x[i];
		y1 = y[i] + xx * (q[i] + xx * (r[i]  + s[i] * xx));

		return y1;
	}
}

  3. JavaScript

    　　ここをクリック（オブジェクトを使用した場合）すると，任意のデータに対して，スプライン補間法による計算結果を画面上で求めることができます． 
    <!DOCTYPE HTML>

<HTML>

<HEAD>

	<TITLE>補間法（３次スプライン関数）</TITLE>
	<META HTTP-EQUIV="Content-Type" CONTENT="text/html; charset=utf-8">
	<SCRIPT TYPE="text/javascript">
		function main()
		{
					// データの設定
			let n  = parseInt(document.getElementById("order").value);
			let sp = parseFloat(document.getElementById("step").value);
			let h  = new Array();
			let b  = new Array();
			let d  = new Array();
			let g  = new Array();
			let u  = new Array();
			let r  = new Array();
			let x  = new Array();
			let y  = new Array();
			let s = (document.getElementById("data").value).split(/ {1,}|\n{1,}/);
			let k = 0;
			for (let i1 = 0; i1 < 2*(n+1); i1 += 2) {
				x[k] = parseFloat(s[i1]);
				y[k] = parseFloat(s[i1+1]);
				k++;
			}
					// 出力
			let str = "";
			let x1 = x[0];
			while (x1 < x[n]+0.01*sp) {
				let y1 = spline(n, x, y, x1, h, b, d, g, u, r);
				str = str + "x " + x1 + " y " + y1 + "\n";
				x1 += sp;
			}
			document.getElementById("ans").value = str;
		}

		/**********************************/
		/* ３次スプライン関数による補間   */
		/*      n : 区間の数              */
		/*      x, y : 点の座標           */
		/*      x1 : 補間値を求める値     */
		/*      h, b, d, g, u, r : 作業域 */
		/*      return : 補間値           */
		/*      coded by Y.Suganuma       */
		/**********************************/
		function spline(n, x, y, x1, h, b, d, g, u, r)
		{
			let i = -1;
			let i1;
			let k;
			let y1;
			let qi;
			let si;
			let xx;
						// 区間の決定
			for (i1 = 1; i1 < n && i < 0; i1++) {
				if (x1 < x[i1])
					i = i1 - 1;
			}
			if (i < 0)
				i = n - 1;
						// ステップ１
			for (i1 = 0; i1 < n; i1++)
				h[i1] = x[i1+1] - x[i1];
			for (i1 = 1; i1 < n; i1++) {
				b[i1] = 2.0 * (h[i1] + h[i1-1]);
				d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1]);
			}
						// ステップ２
			g[1] = h[1] / b[1];
			for (i1 = 2; i1 < n-1; i1++)
				g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
			u[1] = d[1] / b[1];
			for (i1 = 2; i1 < n; i1++)
				u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
						// ステップ３
			k      = (i > 1) ? i : 1;
			r[0]   = 0.0;
			r[n]   = 0.0;
			r[n-1] = u[n-1];
			for (i1 = n-2; i1 >= k; i1--)
				r[i1] = u[i1] - g[i1] * r[i1+1];
						// ステップ４
			xx = x1 - x[i];
			qi = (y[i+1] - y[i]) / h[i] - h[i] * (r[i+1] + 2.0 * r[i]) / 3.0;
			si = (r[i+1] - r[i]) / (3.0 * h[i]);
			y1 = y[i] + xx * (qi + xx * (r[i]  + si * xx));

			return y1;
		}
	</SCRIPT>

</HEAD>

<BODY STYLE="font-size: 130%; background-color: #eeffee;">

	<H2 STYLE="text-align:center"><B>補間法（３次スプライン関数）</B></H2>

	<DL>
		<DT>　　テキストフィールドおよびテキストエリアには，例として，6 点 ( 0, 0 ), ( 1, 1 ), ( 2, 4 ) ( 3, 9 ), ( 4, 16 ), ( 5, 25 ) を使用して，3 次スプライン関数による補間を実行するための値が設定されています．他の問題を実行する場合は，それらを適切に修正してください．
	</DL>

	<DIV STYLE="text-align:center">
		区間の数(n)：<INPUT ID="order" STYLE="font-size: 100%" TYPE="text" SIZE="2" VALUE="5"> 
		ステップ：<INPUT ID="step" STYLE="font-size: 100%;" TYPE="text" SIZE="5" VALUE="0.1"> 
		<BUTTON STYLE="font-size: 100%; background-color: pink" onClick="main()">OK</BUTTON><BR><BR>
		データ（x y，(n+1)個）：<TEXTAREA ID="data" COLS="30" ROWS="15" STYLE="font-size: 100%">0.0 0.0
1.0 1.0
2.0 4.0
3.0 9.0
4.0 16.0
5.0 25.0</TEXTAREA><BR><BR>
		<TEXTAREA ID="ans" COLS="50" ROWS="15" STYLE="font-size: 100%"></TEXTAREA>
	</DIV>

</BODY>

</HTML>
    オブジェクトを使用した場合
    <!DOCTYPE HTML>

<HTML>

<HEAD>

	<TITLE>補間法（３次スプライン関数）</TITLE>
	<META HTTP-EQUIV="Content-Type" CONTENT="text/html; charset=utf-8">
	<SCRIPT TYPE="text/javascript">
		spline = null;

		function main()
		{
					// データの設定
			let n  = parseInt(document.getElementById("order").value);
			let sp = parseFloat(document.getElementById("step").value);
			let x  = new Array();
			let y  = new Array();
			let s = (document.getElementById("data").value).split(/ {1,}|\n{1,}/);
			let k = 0;
			for (let i1 = 0; i1 < 2*(n+1); i1 += 2) {
				x[k] = parseFloat(s[i1]);
				y[k] = parseFloat(s[i1+1]);
				k++;
			}

			spline = new Spline(n, x, y);
					// 出力
			let str = "";
			let x1 = x[0];
			while (x1 < x[n]+0.01*sp) {
				let y1 = spline.value(x1);
				str = str + "x " + x1 + " y " + y1 + "\n";
				x1 += sp;
			}
			document.getElementById("ans").value = str;
		}

		/********************************/
		/* ３次スプライン関数による補間 */
		/*      n1 : 区間の数           */
		/*      x1, y1 : 点の座標       */
		/*      coded by Y.Suganuma     */
		/********************************/
		function Spline(n1, x1, y1)
		{
						// 領域の確保
			let h = new Array();
			let b = new Array();
			let d = new Array();
			let g = new Array();
			let u = new Array();
			this.n = n1;
			this.q = new Array();
			this.s = new Array();
			this.r = new Array();
			this.x = new Array();
			this.y = new Array();
			for (let i1 = 0; i1 <= this.n; i1++) {
				this.x[i1] = x1[i1];
				this.y[i1] = y1[i1];
			}
						// ステップ１
			for (let i1 = 0; i1 < this.n; i1++)
				h[i1] = this.x[i1+1] - this.x[i1];
			for (let i1 = 1; i1 < this.n; i1++) {
				b[i1] = 2.0 * (h[i1] + h[i1-1]);
				d[i1] = 3.0 * ((this.y[i1+1] - this.y[i1]) / h[i1] - (this.y[i1] - this.y[i1-1]) / h[i1-1]);
			}
						// ステップ２
			g[1] = h[1] / b[1];
			for (let i1 = 2; i1 < this.n-1; i1++)
				g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
			u[1] = d[1] / b[1];
			for (let i1 = 2; i1 < this.n; i1++)
				u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
						// ステップ３
			this.r[0]        = 0.0;
			this.r[this.n]   = 0.0;
			this.r[this.n-1] = u[this.n-1];
			for (let i1 = this.n-2; i1 >= 1; i1--)
				this.r[i1] = u[i1] - g[i1] * this.r[i1+1];
						// ステップ４
			for (let i1 = 0; i1 < this.n; i1++) {
				this.q [i1] = (this.y[i1+1] - this.y[i1]) / h[i1] - h[i1] * (this.r[i1+1] + 2.0 * this.r[i1]) / 3.0;
				this.s[i1] = (this.r[i1+1] - this.r[i1]) / (3.0 * h[i1]);
			}
		}

		/******************************/
		/* 補間値の計算               */
		/*      x1 : 補間値を求める値 */
		/*      return : 補間値       */
		/******************************/
		Spline.prototype.value = function(x1)
		{
						// 区間の決定
			let i = -1;
			for (let i1 = 1; i1 < spline.n && i < 0; i1++) {
				if (x1 < spline.x[i1])
					i = i1 - 1;
			}
			if (i < 0)
				i = spline.n - 1;
						// 計算
			let xx = x1 - spline.x[i];
			let y1 = spline.y[i] + xx * (spline.q[i] + xx * (spline.r[i]  + spline.s[i] * xx));

			return y1;
		}
	</SCRIPT>

</HEAD>

<BODY STYLE="font-size: 130%; background-color: #eeffee;">

	<H2 STYLE="text-align:center"><B>補間法（３次スプライン関数）</B></H2>

	<DL>
		<DT>　　テキストフィールドおよびテキストエリアには，例として，6 点 ( 0, 0 ), ( 1, 1 ), ( 2, 4 ) ( 3, 9 ), ( 4, 16 ), ( 5, 25 ) を使用して，3 次スプライン関数による補間を実行するための値が設定されています．他の問題を実行する場合は，それらを適切に修正してください．
	</DL>

	<DIV STYLE="text-align:center">
		区間の数(n)：<INPUT ID="order" STYLE="font-size: 100%" TYPE="text" SIZE="2" VALUE="5"> 
		ステップ：<INPUT ID="step" STYLE="font-size: 100%;" TYPE="text" SIZE="5" VALUE="0.1"> 
		<BUTTON STYLE="font-size: 100%; background-color: pink" onClick="main()">OK</BUTTON><BR><BR>
		データ（x y，(n+1)個）：<TEXTAREA ID="data" COLS="30" ROWS="15" STYLE="font-size: 100%">0.0 0.0
1.0 1.0
2.0 4.0
3.0 9.0
4.0 16.0
5.0 25.0</TEXTAREA><BR><BR>
		<TEXTAREA ID="ans" COLS="50" ROWS="15" STYLE="font-size: 100%"></TEXTAREA>
	</DIV>

</BODY>

</HTML>

  4. PHP 

    <?php

/********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/

					// 設定
	$n = 5;
	$h = array($n);
	$b = array($n);
	$d = array($n);
	$g = array($n);
	$u = array($n);
	$r = array($n+1);
	$x = array($n+1);
	$y = array($n+1);

	$x[0] = 0.0;
	$x[1] = 1.0;
	$x[2] = 2.0;
	$x[3] = 3.0;
	$x[4] = 4.0;
	$x[5] = 5.0;
	$y[0] = 0.0;
	$y[1] = 1.0;
	$y[2] = 4.0;
	$y[3] = 9.0;
	$y[4] = 16.0;
	$y[5] = 25.0;
					// 実行と出力
	$x1 = 0.0;
	$sp = 0.1;
	while ($x1 <= 5.01) {
		$y1 = spline($n, $x, $y, $x1, $h, $b, $d, $g, $u, $r);
		printf("%f %f\n", $x1, $y1);
		$x1 += $sp;
	}

/**********************************/
/* ３次スプライン関数による補間   */
/*      n : 区間の数              */
/*      x, y : 点の座標           */
/*      x1 : 補間値を求める値     */
/*      h, b, d, g, u, r : 作業域 */
/*      return : 補間値           */
/*      coded by Y.Suganuma       */
/**********************************/
function spline($n, $x, $y, $x1, $h, $b, $d, $g, $u, $r)
{
	$i = -1;
					// 区間の決定
	for ($i1 = 1; $i1 < $n && $i < 0; $i1++) {
		if ($x1 < $x[$i1])
			$i = $i1 - 1;
	}
	if ($i < 0)
		$i = $n - 1;
					// ステップ１
	for ($i1 = 0; $i1 < $n; $i1++)
		$h[$i1] = $x[$i1+1] - $x[$i1];
	for ($i1 = 1; $i1 < $n; $i1++) {
		$b[$i1] = 2.0 * ($h[$i1] + $h[$i1-1]);
		$d[$i1] = 3.0 * (($y[$i1+1] - $y[$i1]) / $h[$i1] - ($y[$i1] - $y[$i1-1]) / $h[$i1-1]);
	}
					// ステップ２
	$g[1] = $h[1] / $b[1];
	for ($i1 = 2; $i1 < $n-1; $i1++)
		$g[$i1] = $h[$i1] / ($b[$i1] - $h[$i1-1] * $g[$i1-1]);
	$u[1] = $d[1] / $b[1];
	for ($i1 = 2; $i1 < $n; $i1++)
		$u[$i1] = ($d[$i1] - $h[$i1-1] * $u[$i1-1]) / ($b[$i1] - $h[$i1-1] * $g[$i1-1]);
					// ステップ３
	$k       = ($i > 1) ? $i : 1;
	$r[0]    = 0.0;
	$r[$n]   = 0.0;
	$r[$n-1] = $u[$n-1];
	for ($i1 = $n-2; $i1 >= $k; $i1--)
		$r[$i1] = $u[$i1] - $g[$i1] * $r[$i1+1];
					// ステップ４
	$xx = $x1 - $x[$i];
	$qi = ($y[$i+1] - $y[$i]) / $h[$i] - $h[$i] * ($r[$i+1] + 2.0 * $r[$i]) / 3.0;
	$si = ($r[$i+1] - $r[$i]) / (3.0 * $h[$i]);
	$y1 = $y[$i] + $xx * ($qi + $xx * ($r[$i]  + $si * $xx));

	return $y1;
}

?>
    クラスを使用した場合
    <?php

/********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/

/****************/
/* クラスの定義 */
/****************/
class Spline
{
	private $r = array();
	private $q = array();
	private $s = array();
	private $x = array();
	private $y = array();
	private $n;

	/**************************/
	/* コンストラクタ         */
	/*      n1 : 区間の数     */
	/*      x1, y1 : 点の座標 */
	/**************************/
	function Spline($n1, $x1, $y1)
	{
					// 領域の確保
		$this->n = $n1;

		$h = array();
		$b = array();
		$d = array();
		$g = array();
		$u = array();

		for ($i1 = 0; $i1 <= $this->n; $i1++) {
			$this->x[$i1] = $x1[$i1];
			$this->y[$i1] = $y1[$i1];
		}
					// ステップ１
		for ($i1 = 0; $i1 < $this->n; $i1++)
			$h[$i1] = $this->x[$i1+1] - $this->x[$i1];
		for ($i1 = 1; $i1 < $this->n; $i1++) {
			$b[$i1] = 2.0 * ($h[$i1] + $h[$i1-1]);
			$d[$i1] = 3.0 * (($this->y[$i1+1] - $this->y[$i1]) / $h[$i1] - ($this->y[$i1] - $this->y[$i1-1]) / $h[$i1-1]);
		}
					// ステップ２
		$g[1] = $h[1] / $b[1];
		for ($i1 = 2; $i1 < $this->n-1; $i1++)
			$g[$i1] = $h[$i1] / ($b[$i1] - $h[$i1-1] * $g[$i1-1]);
		$u[1] = $d[1] / $b[1];
		for ($i1 = 2; $i1 < $this->n; $i1++)
			$u[$i1] = ($d[$i1] - $h[$i1-1] * $u[$i1-1]) / ($b[$i1] - $h[$i1-1] * $g[$i1-1]);
					// ステップ３
		$this->r[0]          = 0.0;
		$this->r[$this->n]   = 0.0;
		$this->r[$this->n-1] = $u[$this->n-1];
		for ($i1 = $this->n-2; $i1 >= 1; $i1--)
			$this->r[$i1] = $u[$i1] - $g[$i1] * $this->r[$i1+1];
					// ステップ４
		for ($i1 = 0; $i1 < $this->n; $i1++) {
			$this->q[$i1] = ($this->y[$i1+1] - $this->y[$i1]) / $h[$i1] - $h[$i1] * ($this->r[$i1+1] + 2.0 * $this->r[$i1]) / 3.0;
			$this->s[$i1] = ($this->r[$i1+1] - $this->r[$i1]) / (3.0 * $h[$i1]);
		}
	}

	/******************************/
	/* 補間値の計算               */
	/*      x1 : 補間値を求める値 */
	/*      return : 補間値       */
	/******************************/
	function value($x1)
	{
		$i = -1;
					// 区間の決定
		for ($i1 = 1; $i1 < $this->n && $i < 0; $i1++) {
			if ($x1 < $this->x[$i1])
				$i = $i1 - 1;
		}
		if ($i < 0)
			$i = $this->n - 1;
					// 計算
		$xx = $x1 - $this->x[$i];
		$y1 = $this->y[$i] + $xx * ($this->q[$i] + $xx * ($this->r[$i]  + $this->s[$i] * $xx));

		return $y1;
	}
}

					// 設定
	$n = 5;
	$x = array($n+1);
	$y = array($n+1);

	$x[0] = 0.0;
	$x[1] = 1.0;
	$x[2] = 2.0;
	$x[3] = 3.0;
	$x[4] = 4.0;
	$x[5] = 5.0;
	$y[0] = 0.0;
	$y[1] = 1.0;
	$y[2] = 4.0;
	$y[3] = 9.0;
	$y[4] = 16.0;
	$y[5] = 25.0;

	$spline = new Spline($n, $x, $y);
					// 実行と出力
	$x1 = 0.0;
	$sp = 0.1;
	while ($x1 <= 5.01) {
		$y1 = $spline->value($x1);
		printf("%f %f\n", $x1, $y1);
		$x1 += $sp;
	}

?>

  5. Ruby 

    #*******************************/
# 補間法（３次スプライン関数） */
#      coded by Y.Suganuma     */
#*******************************/

#*********************************/
# ３次スプライン関数による補間   */
#      n : 区間の数              */
#      x, y : 点の座標           */
#      x1 : 補間値を求める値     */
#      h, b, d, g, u, r : 作業域 */
#      return : 補間値           */
#      coded by Y.Suganuma       */
#*********************************/
def spline(n, x, y, x1, h, b, d, g, u, r)
					# 区間の決定
	i = -1
	for i1 in 1 ... n
		if x1 < x[i1]
			i = i1 - 1
		end
		if i >= 0
			break
		end
	end
	if i < 0
		i = n - 1
	end
					# ステップ１
	for i1 in 0 ... n
		h[i1] = x[i1+1] - x[i1]
	end
	for i1 in 1 ... n
		b[i1] = 2.0 * (h[i1] + h[i1-1])
		d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1])
	end
					# ステップ２
	g[1] = h[1] / b[1]
	for i1 in 2 ... n-1
		g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1])
	end
	u[1] = d[1] / b[1]
	for i1 in 2 ... n
		u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1])
	end
					# ステップ３
	k      = (i > 1) ? i : 1
	r[0]   = 0.0
	r[n]   = 0.0
	r[n-1] = u[n-1]
	i1     = n - 2
	while i1 >= k
		r[i1] = u[i1] - g[i1] * r[i1+1]
		i1 -= 1;
	end
					# ステップ４
	xx = x1 - x[i]
	qi = (y[i+1] - y[i]) / h[i] - h[i] * (r[i+1] + 2.0 * r[i]) / 3.0
	si = (r[i+1] - r[i]) / (3.0 * h[i])
	y1 = y[i] + xx * (qi + xx * (r[i]  + si * xx))

	return y1
end

					# 設定
n = 5
h = Array.new(n)
b = Array.new(n)
d = Array.new(n)
g = Array.new(n)
u = Array.new(n)
r = Array.new(n+1)
x = Array[0.0, 1.0, 2.0, 3.0, 4.0, 5.0]
y = Array[0.0, 1.0, 4.0, 9.0, 16.0, 25.0]
				# 実行と出力
x1 = 0.0
sp = 0.1
while x1 <= 5.01
	y1 = spline(n, x, y, x1, h, b, d, g, u, r)
	printf("%f %f\n", x1, y1)
	x1 += sp
end
    クラスを使用した場合
    #*******************************/
# 補間法（３次スプライン関数） */
#      coded by Y.Suganuma     */
#*******************************/

#***************/
# クラスの定義 */
#***************/
class Spline

	#*************************/
	# コンストラクタ         */
	#      n1 : 区間の数     */
	#      x1, y1 : 点の座標 */
	#*************************/
	def initialize(n1, x1, y1)
					# 領域の確保
		h = Array.new(n1)
		b = Array.new(n1)
		d = Array.new(n1)
		g = Array.new(n1)
		u = Array.new(n1)
		@_n = n1
		@_q = Array.new(@_n)
		@_s = Array.new(@_n)
		@_r = Array.new(@_n+1)
		@_x = Array.new(@_n+1)
		@_y = Array.new(@_n+1)

		for i1 in 0 ... @_n+1
			@_x[i1] = x1[i1]
			@_y[i1] = y1[i1]
		end
					# ステップ１
		for i1 in 0 ... @_n
			h[i1] = @_x[i1+1] - @_x[i1]
		end
		for i1 in 1 ... @_n
			b[i1] = 2.0 * (h[i1] + h[i1-1])
			d[i1] = 3.0 * ((@_y[i1+1] - @_y[i1]) / h[i1] - (@_y[i1] - @_y[i1-1]) / h[i1-1])
		end
					# ステップ２
		g[1] = h[1] / b[1]
		for i1 in 2 ... @_n-1
			g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1])
		end
		u[1] = d[1] / b[1]
		for i1 in 2 ... @_n
			u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1])
		end
					# ステップ３
		@_r[0]     = 0.0
		@_r[@_n]   = 0.0
		@_r[@_n-1] = u[@_n-1]
		i1         = @_n - 2
		while i1 >= 1
			@_r[i1] = u[i1] - g[i1] * @_r[i1+1]
			i1     -= 1
		end
					# ステップ４
		for i1 in 0 ... @_n
			@_q[i1] = (@_y[i1+1] - @_y[i1]) / h[i1] - h[i1] * (@_r[i1+1] + 2.0 * @_r[i1]) / 3.0
			@_s[i1] = (@_r[i1+1] - @_r[i1]) / (3.0 * h[i1])
		end
	end

	#*****************************/
	# 補間値の計算               */
	#      x1 : 補間値を求める値 */
	#      return : 補間値       */
	#*****************************/
	def value(x1)
					# 区間の決定
		i = -1
		for i1 in 1 ... @_n
			if (x1 < @_x[i1])
				i = i1 - 1
			end
			if i >= 0
				break
			end
		end
		if (i < 0)
			i = @_n - 1
		end
					# 計算
		xx = x1 - @_x[i]
		y1 = @_y[i] + xx * (@_q[i] + xx * (@_r[i]  + @_s[i] * xx))

		return y1
	end
end

					# 設定
n = 5
x = Array[0.0, 1.0, 2.0, 3.0, 4.0, 5.0]
y = Array[0.0, 1.0, 4.0, 9.0, 16.0, 25.0]

spline = Spline.new(n, x, y)
				# 実行と出力
x1 = 0.0
sp = 0.1
while x1 <= 5.01
	y1 = spline.value(x1)
	printf("%f %f\n", x1, y1)
	x1 += sp
end

  6. Python 

    # -*- coding: UTF-8 -*-
import numpy as np
from math import *

#*********************************/
# ３次スプライン関数による補間   */
#      n : 区間の数              */
#      x, y : 点の座標           */
#      x1 : 補間値を求める値     */
#      h, b, d, g, u, r : 作業域 */
#      return : 補間値           */
#      coded by Y.Suganuma       */
#*********************************/
def spline(n, x, y, x1, h, b, d, g, u, r) :
					# 区間の決定
	i = -1
	for i1 in range(1, n) :
		if x1 < x[i1] :
			i = i1 - 1
		if i >= 0 :
			break
	if i < 0 :
		i = n - 1
					# ステップ１
	for i1 in range(0, n) :
		h[i1] = x[i1+1] - x[i1]
	for i1 in range(1, n) :
		b[i1] = 2.0 * (h[i1] + h[i1-1])
		d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1])
					# ステップ２
	g[1] = h[1] / b[1]
	for i1 in range(2, n-1) :
		g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1])
	u[1] = d[1] / b[1]
	for i1 in range(2, n) :
		u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1])
					# ステップ３
	if i1 > 1 :
		k = i
	else :
		k = 1
	r[0]   = 0.0
	r[n]   = 0.0
	r[n-1] = u[n-1]
	for i1 in range(n-2, k-1, -1) :
		r[i1] = u[i1] - g[i1] * r[i1+1]
					# ステップ４
	xx = x1 - x[i]
	qi = (y[i+1] - y[i]) / h[i] - h[i] * (r[i+1] + 2.0 * r[i]) / 3.0
	si = (r[i+1] - r[i]) / (3.0 * h[i])
	y1 = y[i] + xx * (qi + xx * (r[i]  + si * xx))

	return y1

#*******************************/
# 補間法（３次スプライン関数） */
#      coded by Y.Suganuma     */
#*******************************/
					# 設定
n = 5
x = np.array([0, 1, 2, 3, 4, 5], np.float)
y = np.array([0, 1, 4, 9, 16, 25], np.float)
h = np.empty(n, np.float)
b = np.empty(n, np.float)
d = np.empty(n, np.float)
g = np.empty(n, np.float)
u = np.empty(n, np.float)
r = np.empty(n+1, np.float)
					# 実行と出力
x1 = 0.0
sp = 0.1
while x1 <= 5.01 :
	y1 = spline(n, x, y, x1, h, b, d, g, u, r)
	print(x1, y1)
	x1 += sp
    クラスを使用した場合
    # -*- coding: UTF-8 -*-
import numpy as np
from math import *

############################################
# クラスSplineの定義
#      coded by Y.Suganuma
############################################

class Spline :

	#########################
	# コンストラクタ
	#      n1 : 区間の数
	#      x1, y1 : 点の座標
	#########################

	def __init__(self, n1, x1, y1) :
			# 設定
		h = np.empty(n1, np.float)
		b = np.empty(n1, np.float)
		d = np.empty(n1, np.float)
		g = np.empty(n1, np.float)
		u = np.empty(n1, np.float)
		self.n = n1
		self.q = np.empty(self.n, np.float)
		self.s = np.empty(self.n, np.float)
		self.r = np.empty(self.n+1, np.float)
		self.x = np.empty(self.n+1, np.float)
		self.y = np.empty(self.n+1, np.float)

		for i1 in range(0, self.n+1) :
			self.x[i1] = x1[i1]
			self.y[i1] = y1[i1]
			# ステップ１
		for i1 in range(0, self.n) :
			h[i1] = self.x[i1+1] - self.x[i1]
		for i1 in range(1, self.n) :
			b[i1] = 2.0 * (h[i1] + h[i1-1])
			d[i1] = 3.0 * ((self.y[i1+1] - self.y[i1]) / h[i1] - (self.y[i1] - self.y[i1-1]) / h[i1-1])
			# ステップ２
		g[1] = h[1] / b[1]
		for i1 in range(2, self.n-1) :
			g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1])
		u[1] = d[1] / b[1]
		for i1 in range(2, self.n) :
			u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1])
			# ステップ３
		self.r[0]        = 0.0
		self.r[self.n]   = 0.0
		self.r[self.n-1] = u[self.n-1]
		for i1 in range(self.n-2, 0, -1) :
			self.r[i1] = u[i1] - g[i1] * self.r[i1+1]
			# ステップ４
		for i1 in range(0, self.n) :
			self.q[i1] = (self.y[i1+1] - self.y[i1]) / h[i1] - h[i1] * (self.r[i1+1] + 2.0 * self.r[i1]) / 3.0
			self.s[i1] = (self.r[i1+1] - self.r[i1]) / (3.0 * h[i1])

	#########################
	# 補間値の計算
	#      x1 : 補間値を求める値
	#      return : 補間値
	#########################

	def value(self, x1) :
			# 区間の決定
		i = -1
		for i1 in range(1, self.n) :
			if x1 < self.x[i1] :
				i = i1 - 1
				break
		if i < 0 :
			i = self.n - 1
			# 計算
		xx = x1 - self.x[i]
		y1 = self.y[i] + xx * (self.q[i] + xx * (self.r[i]  + self.s[i] * xx))

		return y1

############################################
# 補間法（３次スプライン関数）
#      coded by Y.Suganuma
############################################

			# 設定
n = 5
x = np.array([0, 1, 2, 3, 4, 5], np.float)
y = np.array([0, 1, 4, 9, 16, 25], np.float)

spline = Spline(n, x, y)
			# 実行と出力
x1 = 0.0
sp = 0.1
while x1 <= 5.01 :
	y1 = spline.value(x1)
	print(x1, y1)
	x1 += sp

  7. C# 

    /********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
using System;

class Program
{
	static void Main()
	{
		Test1 ts = new Test1();
	}
}

class Test1
{
	public Test1()
	{
					// 設定
		int n = 5;
		double[] h = new double [n];
		double[] b = new double [n];
		double[] d = new double [n];
		double[] g = new double [n];
		double[] u = new double [n];
		double[] r = new double [n+1];
		double[] x = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0};
		double[] y = {0.0, 1.0, 4.0, 9.0, 16.0, 25.0};
					// 実行と出力
		double x1 = 0.0;
		double sp = 0.1;
		while (x1 <= 5.01) {
			double y1 = spline(n, x, y, x1, h, b, d, g, u, r);
			Console.WriteLine(x1 + " " + y1);
			x1 += sp;
		}
	}

	/**********************************/
	/* ３次スプライン関数による補間   */
	/*      n : 区間の数              */
	/*      x, y : 点の座標           */
	/*      x1 : 補間値を求める値     */
	/*      h, b, d, g, u, r : 作業域 */
	/*      return : 補間値           */
	/*      coded by Y.Suganuma       */
	/**********************************/
	double spline(int n, double[] x, double[] y, double x1, double[] h, double[] b,
	              double[] d, double[] g, double[] u, double[] r)
	{
					// 区間の決定
		int i = -1;
		for (int i1 = 1; i1 < n && i < 0; i1++) {
			if (x1 < x[i1])
				i = i1 - 1;
		}
		if (i < 0)
			i = n - 1;
					// ステップ１
		for (int i1 = 0; i1 < n; i1++)
			h[i1] = x[i1+1] - x[i1];
		for (int i1 = 1; i1 < n; i1++) {
			b[i1] = 2.0 * (h[i1] + h[i1-1]);
			d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1]);
		}
					// ステップ２
		g[1] = h[1] / b[1];
		for (int i1 = 2; i1 < n-1; i1++)
			g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
		u[1] = d[1] / b[1];
		for (int i1 = 2; i1 < n; i1++)
			u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
					// ステップ３
		int k  = (i > 1) ? i : 1;
		r[0]   = 0.0;
		r[n]   = 0.0;
		r[n-1] = u[n-1];
		for (int i1 = n-2; i1 >= k; i1--)
			r[i1] = u[i1] - g[i1] * r[i1+1];
					// ステップ４
		double xx = x1 - x[i];
		double qi = (y[i+1] - y[i]) / h[i] - h[i] * (r[i+1] + 2.0 * r[i]) / 3.0;
		double si = (r[i+1] - r[i]) / (3.0 * h[i]);
		double y1 = y[i] + xx * (qi + xx * (r[i]  + si * xx));

		return y1;
	}
}
    クラスを使用した場合
    /********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
using System;

class Program
{
	static void Main()
	{
					// 設定
		int n = 5;
		double[] x = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0};
		double[] y = {0.0, 1.0, 4.0, 9.0, 16.0, 25.0};

		Spline spline = new Spline(n, x, y);
					// 実行と出力
		double x1 = 0.0;
		double sp = 0.1;
		while (x1 <= 5.01) {
			double y1 = spline.value(x1);
			Console.WriteLine(x1 + " " + y1);
			x1 += sp;
		}
	}
}

/********************************/
/* 補間法（３次スプライン関数） */
/*      coded by Y.Suganuma     */
/********************************/
class Spline
{
	private double[] q;
	private double[] s;
	private double[] r;
	private double[] x;
	private double[] y;
	private int n;

	/**************************/
	/* コンストラクタ         */
	/*      n1 : 区間の数     */
	/*      x1, y1 : 点の座標 */
	/**************************/
	public Spline(int n1, double[] x1, double[] y1)
	{
						// 領域の確保
		n = n1;

		double[] h = new double [n];
		double[] b = new double [n];
		double[] d = new double [n];
		double[] g = new double [n];
		double[] u = new double [n];
		q = new double [n];
		s = new double [n];
		r = new double [n+1];
		x = new double [n+1];
		y = new double [n+1];

		for (int i1 = 0; i1 <= n; i1++) {
			x[i1] = x1[i1];
			y[i1] = y1[i1];
		}
						// ステップ１
		for (int i1 = 0; i1 < n; i1++)
			h[i1] = x[i1+1] - x[i1];
		for (int i1 = 1; i1 < n; i1++) {
			b[i1] = 2.0 * (h[i1] + h[i1-1]);
			d[i1] = 3.0 * ((y[i1+1] - y[i1]) / h[i1] - (y[i1] - y[i1-1]) / h[i1-1]);
		}
						// ステップ２
		g[1] = h[1] / b[1];
		for (int i1 = 2; i1 < n-1; i1++)
			g[i1] = h[i1] / (b[i1] - h[i1-1] * g[i1-1]);
		u[1] = d[1] / b[1];
		for (int i1 = 2; i1 < n; i1++)
			u[i1] = (d[i1] - h[i1-1] * u[i1-1]) / (b[i1] - h[i1-1] * g[i1-1]);
						// ステップ３
		r[0]   = 0.0;
		r[n]   = 0.0;
		r[n-1] = u[n-1];
		for (int i1 = n-2; i1 >= 1; i1--)
			r[i1] = u[i1] - g[i1] * r[i1+1];
						// ステップ４
		for (int i1 = 0; i1 < n; i1++) {
			q[i1] = (y[i1+1] - y[i1]) / h[i1] - h[i1] * (r[i1+1] + 2.0 * r[i1]) / 3.0;
			s[i1] = (r[i1+1] - r[i1]) / (3.0 * h[i1]);
		}
	}

	/******************************/
	/* 補間値の計算               */
	/*      x1 : 補間値を求める値 */
	/*      return : 補間値       */
	/******************************/
	public double value(double x1)
	{
						// 区間の決定
		int i = -1;
		for (int i1 = 1; i1 < n && i < 0; i1++) {
			if (x1 < x[i1])
				i = i1 - 1;
		}
		if (i < 0)
			i = n - 1;
						// 計算
		double xx = x1 - x[i];
		double y1 = y[i] + xx * (q[i] + xx * (r[i]  + s[i] * xx));

		return y1;
	}
}

  8. VB 

    '******************************'
' 補間法（３次スプライン関数） '
'      coded by Y.Suganuma     '
'******************************'
Module Test

	Sub Main()
					' 設定
		Dim n As Integer = 5
		Dim h(n) As Double
		Dim b(n) As Double
		Dim d(n) As Double
		Dim g(n) As Double
		Dim u(n) As Double
		Dim r(n+1) As Double
		Dim x() As Double = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0}
		Dim y() As Double = {0.0, 1.0, 4.0, 9.0, 16.0, 25.0}
					' 実行と出力
		Dim x1 As Double = 0.0
		Dim sp As Double = 0.1
		Do While x1 <= 5.01
			Dim y1 As Double = spline(n, x, y, x1, h, b, d, g, u, r)
			Console.WriteLine(x1 & " " & y1)
			x1 += sp
		Loop

	End Sub

	'********************************'
	' ３次スプライン関数による補間   '
	'      n : 区間の数              '
	'      x, y : 点の座標           '
	'      x1 : 補間値を求める値     '
	'      h, b, d, g, u, r : 作業域 '
	'      return : 補間値           '
	'      coded by Y.Suganuma       '
	'********************************'
	Function spline(n As Integer, x() As Double, y() As Double, x1 As Double,
	                h() As Double, b() As Double, d() As Double, g() As Double,
	                u() As Double, r() As Double)

					' 区間の決定
		Dim i As Integer  = -1
		Dim i1 As Integer = 1
		Do While i1 < n and i < 0
			If x1 < x(i1)
				i = i1 - 1
			End If
			i1 += 1
		Loop
		If i < 0
			i = n - 1
		End If
					' ステップ１
		For i1 = 0 To n-1
			h(i1) = x(i1+1) - x(i1)
		Next
		For i1 = 1 To n-1
			b(i1) = 2.0 * (h(i1) + h(i1-1))
			d(i1) = 3.0 * ((y(i1+1) - y(i1)) / h(i1) - (y(i1) - y(i1-1)) / h(i1-1))
		Next
					' ステップ２
		g(1) = h(1) / b(1)
		For i1 = 2 To n-2
			g(i1) = h(i1) / (b(i1) - h(i1-1) * g(i1-1))
		Next
		u(1) = d(1) / b(1)
		For i1 = 2 To n-1
			u(i1) = (d(i1) - h(i1-1) * u(i1-1)) / (b(i1) - h(i1-1) * g(i1-1))
		Next
					' ステップ３
		Dim k As Integer = i
		If i <= 1
			k = 1
		End If
		r(0)   = 0.0
		r(n)   = 0.0
		r(n-1) = u(n-1)
		For i1 = n-2 To k Step -1
			r(i1) = u(i1) - g(i1) * r(i1+1)
		Next
					' ステップ４
		Dim xx As Double = x1 - x(i)
		Dim qi As Double = (y(i+1) - y(i)) / h(i) - h(i) * (r(i+1) + 2.0 * r(i)) / 3.0
		Dim si As Double = (r(i+1) - r(i)) / (3.0 * h(i))
		Dim y1 As Double = y(i) + xx * (qi + xx * (r(i)  + si * xx))

		Return y1

	End Function

End Module
    クラスを使用した場合
    '******************************'
' 補間法（３次スプライン関数） '
'      coded by Y.Suganuma     '
'******************************'
Module Test

	Sub Main()
					' 設定
		Dim n As Integer  = 5
		Dim x() As Double = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0}
		Dim y() As Double = {0.0, 1.0, 4.0, 9.0, 16.0, 25.0}
	
		Dim sp As Spline = new Spline(n, x, y)
					' 実行と出力
		Dim x1 As Double  = 0.0
		Dim sp1 As Double = 0.1
		Do While x1 <= 5.01
			Dim y1 As Double = sp.value(x1)
			Console.WriteLine(x1 & " " & y1)
			x1 += sp1
		Loop
	End Sub

	'******************************'
	' 補間法（３次スプライン関数） '
	'      coded by Y.Suganuma     '
	'******************************'
	Class Spline

		Private n As Integer
		Private q() As Double
		Private s() As Double
		Private r() As Double
		Private x() As Double
		Private y() As Double
	
		'************************'
		' コンストラクタ         '
		'      n1 : 区間の数     '
		'      x1, y1 : 点の座標 '
		'************************'
		Public Sub New (n1 As Integer, x1() As Double, y1() As Double)
							' 領域の確保
			n = n1
			ReDim q(n)
			ReDim s(n)
			ReDim r(n+1)
			ReDim x(n+1)
			ReDim y(n+1)
	
			Dim h(n) As Double
			Dim b(n) As Double
			Dim d(n) As Double
			Dim g(n) As Double
			Dim u(n) As Double
	
			For i1 As Integer = 0 To n
				x(i1) = x1(i1)
				y(i1) = y1(i1)
			Next
 							' ステップ１
			For i1 As Integer = 0 To n-1
				h(i1) = x(i1+1) - x(i1)
			Next
			For i1 As Integer = 1 To n-1
				b(i1) = 2.0 * (h(i1) + h(i1-1))
				d(i1) = 3.0 * ((y(i1+1) - y(i1)) / h(i1) - (y(i1) - y(i1-1)) / h(i1-1))
			Next
							' ステップ２
			g(1) = h(1) / b(1)
			For i1 As Integer = 2 To n-2
				g(i1) = h(i1) / (b(i1) - h(i1-1) * g(i1-1))
			Next
			u(1) = d(1) / b(1)
			For i1 As Integer = 2 To n-1
				u(i1) = (d(i1) - h(i1-1) * u(i1-1)) / (b(i1) - h(i1-1) * g(i1-1))
			Next
							' ステップ３
			r(0)   = 0.0
			r(n)   = 0.0
			r(n-1) = u(n-1)
			For i1 As Integer = n-2 To 1 Step -1
				r(i1) = u(i1) - g(i1) * r(i1+1)
			Next
							' ステップ４
			For i1 As Integer = 0 To n-1
				q(i1) = (y(i1+1) - y(i1)) / h(i1) - h(i1) * (r(i1+1) + 2.0 * r(i1)) / 3.0
				s(i1) = (r(i1+1) - r(i1)) / (3.0 * h(i1))
			Next

		End Sub
	
		'****************************'
		' 補間値の計算               '
		'      x1 : 補間値を求める値 '
		'      return : 補間値       '
		'****************************'
		Function value(x1 As Double) As Double
							' 区間の決定
			Dim i As Integer  = -1
			Dim i1 As Integer = 1
			Do While i1 < n and i < 0
				If x1 < x(i1)
					i = i1 - 1
				End If
				I1 += 1
			Loop
			If i < 0
				i = n - 1
			End If
							' 計算
			Dim xx As Double = x1 - x(i)
			Dim y1 As Double = y(i) + xx * (q(i) + xx * (r(i)  + s(i) * xx))
	
			Return y1
	
		End Function
	
	End Class

End Module

情報学部 菅沼ホーム 目次 索引