Levenberg-Marquardt c++程序

Levenberg-Marquardt算法是一种用于非线性最小二乘问题的优化算法。它结合了梯度下降法和牛顿法的优点，通过迭代更新参数来最小化目标函数。在C++中实现Levenberg-Marquardt算法通常涉及以下几个步骤：

1. 定义目标函数和其雅可比矩阵（Jacobian matrix）。
2. 初始化参数向量和阻尼因子。
3. 计算目标函数的梯度和雅可比矩阵。
4. 更新参数向量，使用阻尼因子调整步长。
5. 检查收敛条件，如果未满足则重复步骤3-4。
6. 输出最终的参数估计值。

下面是一个简单的C++程序示例，演示了如何使用Levenberg-Marquardt算法来拟合一个二次函数：

#include <iostream>
#include <Eigen/Dense>

using namespace Eigen;
using namespace std;

// 目标函数：f(x) = a*x^2 + b*x + c
VectorXd func(const VectorXd &params, const VectorXd &x) {
    return params[0] * x.array().square() + params[1] * x.array() + params[2];
}

// 目标函数的梯度
VectorXd gradient(const VectorXd &params, const VectorXd &x, const VectorXd &y) {
    VectorXd residuals = func(params, x) - y;
    VectorXd grad = VectorXd::Zero(params.size());
    grad[0] = (residuals.array() * x.array().square()).sum();
    grad[1] = (residuals.array() * x.array()).sum();
    grad[2] = residuals.sum();
    return grad;
}

// 目标函数的雅可比矩阵
MatrixXd jacobian(const VectorXd &params, const VectorXd &x) {
    MatrixXd jac = MatrixXd::Zero(x.size(), params.size());
    for (int i = 0; i < x.size(); ++i) {
        jac(i, 0) = x[i] * x[i];
        jac(i, 1) = x[i];
        jac(i, 2) = 1;
    }
    return jac;
}

// Levenberg-Marquardt算法
VectorXd levenbergMarquardt(const VectorXd &initial_params, const VectorXd &x, const VectorXd &y, double lambda) {
    VectorXd params = initial_params;
    for (int iter = 0; iter < 100; ++iter) {
        MatrixXd J = jacobian(params, x);
        VectorXd g = gradient(params, x, y);
        MatrixXd H = J.transpose() * J;
        VectorXd h = J.transpose() * g;
        H += lambda * MatrixXd::Identity(params.size(), params.size());
        VectorXd delta = H.ldlt().solve(h);
        params -= delta;
        if (delta.norm() < 1e-6) break;
    }
    return params;
}

int main() {
    VectorXd x(10), y(10);
    for (int i = 0; i < 10; ++i) {
        x[i] = i;
        y[i] = 1.0 * i * i + 2.0 * i + 3.0; // y = x^2 + 2x + 3
    }
    VectorXd initial_params(3);
    initial_params << 0, 0, 0; // 初始猜测参数
    double lambda = 0.01; // 阻尼因子
    VectorXd estimated_params = levenbergMarquardt(initial_params, x, y, lambda);
    cout << "Estimated parameters: " << estimated_params.transpose() << endl;
    return 0;
}

这个程序首先定义了一个二次函数及其梯度和雅可比矩阵，然后实现了Levenberg-Marqu