補間法(ラグランジュ)

/****************************/
/* 補間法(ラグランジュ)   */
/*      coded by Y.Suganuma */
/****************************/
import java.io.*;

public class Test {
	public static void main(String args[]) throws IOException
	{
		double x[], x1, y[], y1;
		int i1, n;

		x = new double [3];
		y = new double [3];
					// 線形補間
		n    = 1;
		x[0] = 0.0;
		x[1] = 2.0;
		y[0] = 0.0;
		y[1] = 4.0;

		System.out.println("線形補間");
		x1 = 0.0;
		for (i1 = 0; i1 <= 8; i1++) {
			y1 = lagrange(n, x, y, x1);
			System.out.println("   x " + x1 + " y " + y1);
			x1 += 0.25;
		}

					// 2次補間
		n    = 2;
		x[0] = 0.0;
		x[1] = 1.0;
		x[2] = 2.0;
		y[0] = 0.0;
		y[1] = 1.0;
		y[2] = 4.0;

		System.out.println("2次補間\n");
		x1 = 0.0;
		for (i1 = 0; i1 <= 8; i1++) {
			y1 = lagrange(n, x, y, x1);
			System.out.println("   x " + x1 + " y " + y1);
			x1 += 0.25;
		}
	}

	/**********************************/
	/* ラグランジュ補間法             */
	/*      n : 次数(n=1は線形補間) */
	/*      x, y : (n+1)個の点の座標  */
	/*      x1 : 補間したい点のx座標 */
	/*      return : 補間結果         */
	/**********************************/
	static double lagrange(int n, double x[], double y[], double x1)
	{
		double s1, s2, y1 = 0.0;
		int i1, i2;

		for (i1 = 0; i1 <= n; i1++) {
			s1 = 1.0;
			s2 = 1.0;
			for (i2 = 0; i2 <= n; i2++) {
				if (i2 != i1) {
					s1 *= (x1 - x[i2]);
					s2 *= (x[i1] - x[i2]);
				}
			}
			y1 += y[i1] * s1 / s2;
		}

		return y1;
	}
}