拉格朗日插值算法

public class Lagrange
{
    /// <summary>
    /// X各点坐标组成的数组
    /// </summary>
    public double[] x { get; set; }

    /// <summary>
    /// X各点对应的Y坐标值组成的数组
    /// </summary>
    public double[] y { get; set; }

    /// <summary>
    /// x数组或者y数组中元素的个数, 注意两个数组中的元素个数需要一样
    /// </summary>
    public int itemNum { get; set; }

    /// <summary>
    /// 初始化拉格朗日插值
    /// </summary>
    /// <param name="x">X各点坐标组成的数组</param>
    /// <param name="y">X各点对应的Y坐标值组成的数组</param>
    public Lagrange(double[] x, double[] y)
    {
        this.x = x; this.y = y;
        this.itemNum = x.Length;
    }

    /// <summary>
    /// 获得某个横坐标对应的Y坐标值
    /// </summary>
    /// <param name="xValue">x坐标值</param>
    /// <returns></returns>
    public double GetValue(double xValue)
    {
        //用于累乘数组始末下标
        int start, end;
        //返回值
        double value = 0.0;
        //如果初始的离散点为空, 返回0
        if (itemNum < 1) { return value; }
        //如果初始的离散点只有1个, 返回该点对应的Y值
        if (itemNum == 1) { value = y[0]; return value; }
        //如果初始的离散点只有2个, 进行线性插值并返回插值
        if (itemNum == 2)
        {
            value = (y[0] * (xValue - x[1]) - y[1] * (xValue - x[0])) / (x[0] - x[1]);
            return value;
        }
        //如果插值点小于第一个点X坐标, 取数组前3个点做插值
        if (xValue <= x[1]) { start = 0; end = 2; }
        //如果插值点大于等于最后一个点X坐标, 取数组最后3个点做插值
        else if (xValue >= x[itemNum - 2]) { start = itemNum - 3; end = itemNum - 1; }
        //除了上述的一些特殊情况, 通常情况如下
        else
        {
            start = 1; end = itemNum;
            int temp;
            //使用二分法决定选择哪三个点做插值
            while ((end - start) != 1)
            {
                temp = (start + end) / 2;
                if (xValue < x[temp - 1])
                    end = temp;
                else
                    start = temp;
            }
            start--; end--;
            //看插值点跟哪个点比较靠近
            if (Math.Abs(xValue - x[start]) < Math.Abs(xValue - x[end]))
                start--;
            else
                end++;
        }
        //这时已经确定了取哪三个点做插值, 第一个点为x[start]
        double valueTemp;
        //注意是二次的插值公式
        for (int i = start; i <= end; i++)
        {
            valueTemp = 1.0;
            for (int j = start; j <= end; j++)
                if (j != i)
                    valueTemp *= (double)(xValue - x[j]) / (double)(x[i] - x[j]);
            value += valueTemp * y[i];
        }
        return value;
    }
}

调用实例

double[] xValues = new double[] { 0.00, 40.00, 80.00, 120.00, 160.00, 200.00, 240.00, 280.00, 320.00, 360.00, 400.00, 440.00, 480.00, 520.00, 560.00, 600.00, 640.00, 680.00, 720.00, 760.00, 800.00, 840.00, 880.00, 920.00, 960.00, 1000.00, 1040.00, 1080.00, 1120.00, 1160.00, 1200.00, 1240.00, 1280.00, 1320.00, 1360.00, 1400.00, 1440.00, 1480.00, 1520.00, 1560.00, 1600.00, 1640.00, 1680.00, 1720.00, 1760.00, 1800.00, 1840.00, 1880.00, 1920.00, 1960.00, 2000.00, 2040.00, 2080.00, 2120.00, 2160.00, 2200.00, 2240.00, 2280.00, 2320.00, 2360.00, 2400.00, 2440.00, 2480.00, 2520.00, 2560.00, 2600.00, 2640.00, 2680.00, 2720.00, 2760.00, 2800.00 };
double[] yValues = new double[] { 0.00, 148.00, 420.00, 773.00, 1189.00, 1661.00, 2180.00, 2743.00, 3344.00, 3982.00, 4652.00, 5353.00, 6081.00, 6836.00, 7616.00, 8418.00, 9241.00, 10084.00, 10945.00, 11824.00, 12718.00, 13627.00, 14549.00, 15483.00, 16429.00, 17385.00, 18351.00, 19324.00, 20305.00, 21292.00, 22285.00, 23282.00, 24283.00, 25287.00, 26292.00, 27298.00, 28305.00, 29310.00, 30313.00, 31314.00, 32312.00, 33304.00, 34292.00, 35273.00, 36246.00, 37212.00, 38168.00, 39113.00, 40048.00, 40970.00, 41879.00, 42773.00, 43651.00, 44513.00, 45355.00, 46179.00, 46981.00, 47760.00, 48515.00, 49244.00, 49945.00, 50615.00, 51252.00, 51854.00, 52417.00, 52936.00, 53407.00, 53824.00, 54176.00, 54448.00, 54597.00 };


Lagrange lagrange = new Lagrange(xValues, yValues);
double xx = lagrange.GetValue(281.0);