输入:插值节点数组、插值节点处的函数值数组,待求点
输出:函数值
代码如下:把printf的注释取消掉,能打印出中间计算过程,包括Lagrange多项式的求解,多项式每一项等等(代码多次修改,这些prinft不知道还有没有疏漏,日后再检查检查)
注:如果要求解插值区间内的多个点的函数值,把newx改成数组类型,把result改成指针类型,然后加一层for循环即可
/* 这里形参必须是int x[], 写作int[] x产生语法错误 n次插值Lagrange形式 x: 插值节点 y : 插值节点函数值 len : 插值节点个数 newX: 所求节点 */ double lagrangeInterpolation(double x[], double y[], int len, double newx) { //printf("-------------- "); double result = 0;// 这里要记得初始化,否则结果错 double L; // lagrange interpolation polyminoal for (int j = 0; j < len; j++) { // every interpolation polyminoal L = y[j]; for (int k = 0; k < len; k++) { // every term if (j == k) { continue; } //printf("newx[i]=%lf, x[j]=%lf, x[k]=%lf, (newx[i] - x[k]) / (x[j] - x[k])=%lf ", newx[i], x[j], x[k], (newx[i] - x[k]) / (x[j] - x[k])); //printf("L*(newx[i] - x[k]) / (x[j] - x[k])=%lf ", L*(newx[i] - x[k]) / (x[j] - x[k])); L *= (newx - x[k]) / (x[j] - x[k]); } //printf("result[%d]=%lf ", i, result[i]); result += L; //printf("result[%d]=%lf , iter=%d ", i, result[i], j); //printf("array[%d] is: %lf ", i, result[i]); } return result; }
主函数示例:
int main() { printf("Nuerical Analysis Lagrange Interpolation! "); // 插值节点与函数值double x[] = { 0.46, 0.47, 0.48, 0.49 }; int len = sizeof(x) / sizeof(x[0]); double y[] = { 0.4846555, 0.4937452, 0.5027498, 0.5116683 }; // 待求节点 double newX = 0.472; double result = lagrangeInterpolation(x, y, len, newX); printf("f(0.472) = %lf ", result);return 0; }
结果0.495553