4

我正在从事一些科学项目,我需要广义缩减梯度算法的 C 语言实现来进行非线性优化。是否有任何库或只是一段代码?或者,请为非线性多变量问题提出任何其他解决方案。我正在寻找使用 4 个自变量和 2 个常数来构建优化模型:该模型是非线性的。我已经检查了 Microsoft Excel 的 Solver,使用广义缩减梯度 (GRG) 可以完美地解决这个模型,但我需要用 C 语言进行模拟。

这是我的 excel 解决方案: http : //speedy.sh/SEdZj/eof-cs-rest.xlsm 我使用带有 GRG 算法的 Microsoft Excel Solver 来搜索 SS 的最小值,输出是 Const_a 和 Const_b 的值。

4

1 回答 1

1

由GAMS分发的 CONOPT似乎是 GRG 的既定实现,但不是免费的(尽管演示对您来说可能就足够了)。

Alglib 在这里实现了非线性 Levenberg-Marquardt 算法,并且是 GPL / 商业许可的。

下面使用 alglib 的示例代码:

/*
 * Simple optimiser example
 *
 * nl_opt.cpp
 *
 * Compile with eg 'g++ -I../tools/alglib/src ../tools/alglib/src/ap.cpp ../tools/alglib/src/alglibinternal.cpp ../tools/alglib/src/linalg.cpp ../tools/alglib/src/alglibmisc.cpp ../tools/alglib/src/solvers.cpp ../tools/alglib/src/optimization.cpp nl_opt.cpp -o opt'
 *
 */

#include "stdafx.h"
#include <iostream>
#include <cmath>

#include "optimization.h"

using namespace std;

double fn(double a1, double a2, double a3, double x, double A, double B)
{
    return A * exp(-x*(a1*B*B+a2*B+a3));
}

struct problem
{
    double *m_a1s;
    double *m_a2s;
    double *m_a3s;
    double *m_xs;
    double *m_ys;

    int m_n;

    problem(double *a1s, double *a2s, double *a3s, double *xs, double *ys, int n) 
        : m_a1s(a1s), m_a2s(a2s), m_a3s(a3s), m_xs(xs), m_ys(ys), m_n(n)
    {
    }

    void fn_vec(const alglib::real_1d_array &c_var, alglib::real_1d_array &fi, void *ptr)
    {
        double sum = 0.0;
        for(int i = 0; i < m_n; ++i)
        {
            double yhat = fn(m_a1s[i], m_a2s[i], m_a3s[i], m_xs[i], c_var[0], c_var[1]);
            double err_sq = (m_ys[i] - yhat) * (m_ys[i] - yhat);
            sum += err_sq;
        }
        fi[0] = sum;
    }
};

problem *g_p;

void fn_vec(const alglib::real_1d_array &c_var, alglib::real_1d_array &fi, void *ptr)
{
    g_p->fn_vec(c_var, fi, ptr);
}

int main()
{
    cout << "Testing non-linear optimizer..." << endl;

    int n = 5;
    double a1s[] = {2.42, 4.78, 7.25, 9.55, 11.54};
    double a2s[] = {4.25, 5.27, 6.33, 7.32, 8.18};
    double a3s[] = {3.94, 4.05, 4.17, 4.28, 4.37};

    double xs[] = {0.024, 0.036, 0.048, 0.06, 0.072};
    double ys[] = {80, 70, 50, 40, 45};

    double initial[] = {150, 1.75};
    double ss_init = 0.0;

    cout << "Initial problem:" << endl;
    for(int i = 0; i < n; ++i)
    {
        double yhat = fn(a1s[i], a2s[i], a3s[i], xs[i], initial[0], initial[1]);
        double err_sq = (ys[i] - yhat) * (ys[i] - yhat);
        ss_init += err_sq;
        cout << a1s[i] << "\t" << a2s[i] << "\t" << a3s[i] << "\t" 
            << xs[i] << "\t" << ys[i] << "\t" << yhat << "\t" << err_sq << endl;
    }
    cout << "Error: " << ss_init << endl;

    // create problem
    problem p(a1s, a2s, a3s, xs, ys, n);
    g_p = &p;

    // setup solver
    alglib::real_1d_array x = "[150.0, 1.75]";
    double epsg = 0.00000001;
    double epsf = 0;
    double epsx = 0;

    alglib::ae_int_t maxits = 0;
    alglib::minlmstate state;
    alglib::minlmreport report;

    alglib::minlmcreatev(2, x, 0.0001, state);
    alglib::minlmsetcond(state, epsg, epsf, epsx, maxits);

    // optimize
    alglib::minlmoptimize(state, fn_vec);

    alglib::minlmresults(state, x, report);

    cout << "Results:" << endl;

    cout << report.terminationtype << endl;
    cout << x.tostring(2).c_str() << endl;

    double ss_end = 0.0;
    for(int i = 0; i < n; ++i)
    {
        double yhat = fn(a1s[i], a2s[i], a3s[i], xs[i], x[0], x[1]);
        double err_sq = (ys[i] - yhat) * (ys[i] - yhat);
        ss_end += err_sq;
        cout << a1s[i] << "\t" << a2s[i] << "\t" << a3s[i] << "\t"
             << xs[i] << "\t" << ys[i] << "\t" << yhat << "\t" << err_sq << endl;
    }
    cout << "Error: " << ss_end << endl;

    return 0;
}

这给出了示例输出:

./opt 
Testing non-linear optimizer...
Initial problem:
2.42    4.25    3.94    0.024   80  95.5553 241.968
4.78    5.27    4.05    0.036   70  54.9174 227.485
7.25    6.33    4.17    0.048   50  24.8537 632.338
9.55    7.32    4.28    0.06    40  9.3038  942.257
11.54   8.18    4.37    0.072   45  3.06714 1758.36
Error: 3802.41
Results:
2
[92.22,0.57]
2.42    4.25    3.94    0.024   80  77.6579 5.48528
4.78    5.27    4.05    0.036   70  67.599  5.76475
7.25    6.33    4.17    0.048   50  56.6216 43.8456
9.55    7.32    4.28    0.06    40  46.0026 36.0314
11.54   8.18    4.37    0.072   45  36.6279 70.0922
Error: 161.219
于 2013-02-10T14:09:31.840 回答