（十三）计算机数值方法之单纯形搜索法求极值

最新推荐文章于 2024-10-15 19:28:18 发布

热烈又平淡

最新推荐文章于 2024-10-15 19:28:18 发布

阅读量321

点赞数 8

分类专栏：计算机数值方法文章标签：算法 c++

本文链接：https://blog.csdn.net/szt12345__/article/details/139727337

版权

计算机数值方法专栏收录该内容

13 篇文章 0 订阅

订阅专栏

数学问题：

求Rosenbrock函数： $f(x)=\sum_{i=1}^{n-1}\biggl[100(x_{i+1}-x_i^2)^2+(1-x_i)^2\biggr]$ ，当n=5时的最优值，取值范围为 $-30\leq x_i\leq30$ 。初始迭代值为x0=[1,2,3,4,5]，初始步长为h=2，扩张系数alpha=2，收缩系数beta=0.5，压缩系数theta=0.5，误差eps=1e-6。

解决代码：

#include<iostream>  
#include<math.h>  
#include<iomanip>
#include<vector>
using namespace std;
double f(vector<double> x)
{
    double y = 0;
    for (int i = 0; i < 4; i++)
    {
        y = y + 100 * (x[i + 1] - x[i] * x[i]) * (x[i + 1] - x[i] * x[i]) + (1 - x[i]) *(1 - x[i]);
    }
    return y;
}

void SimpleSearch(vector<double>& x0, double h, int n, double alpha, double beta, double theta, double eps, vector<double>& x_opt, double& f_opt, int& stepnum)
{

    vector<vector<double> > xx(n + 1, vector<double>(n));
    for (int j = 0; j < n; j++)
        xx[0][j] = x0[j];
    for (int i = 1; i < n + 1; i++)
    {
        for (int j = 0; j < n; j++)
        {
            xx[i][j] = xx[0][j];
            if (i - 1 == j)
                xx[i][j] = xx[i][j] + h;
        }
    }

    vector<double> Pbest(n);
    vector<double> Pworst(n);
    vector<double> Plworst(n);
    vector<double> Pc(n);
    vector<double> Pr(n);
    vector<double> Pe(n);
    vector<double> Py(n);
    vector<double> Ps(n);
    stepnum = 0;
    value* v = new value[n + 1];
    value vt;
    while (1)
    {
        stepnum = stepnum + 1;

        for (int i = 0; i < n + 1; i++) {
            vt.index = i;
            vt.va = f(xx[i]);
            v[i] = vt;
        }
        double max = 0;
        int ind = 0;
        for (int i = n; i >= 0; i--) {

            max = v[0].va;
            ind = 0;
            for (int j = 0; j <= i; j++) {
                if (v[j].va > max) {
                    max = v[j].va;
                    ind = j;
                }
            }
            vt = v[ind];
            v[ind] = v[i];
            v[i] = vt;
        }
        for (int i = 0; i < n; i++) {

            Pbest[i] = xx[v[0].index][i];
            Pworst[i] = xx[v[n].index][i];
            Plworst[i] = xx[v[n - 1].index][i];
        }

        for (int i = 0; i < n; i++) {
            Pc[i] = 0;
        }
        for (int i = 0; i < n; i++) {

            for (int j = 0; j < n + 1; j++) {
                Pc[i] += xx[j][i];
            }
        }
        for (int i = 0; i < n; i++) {
            Pc[i] -= Pworst[i];
            Pc[i] /= n;
        }
        double error = abs(f(Pworst) - f(Pbest)) / f(Pbest);
        if (error < eps)
        {
            x_opt = Pbest;
            f_opt = f(Pbest);
            break;
        }      
        for (int i = 0; i < n; i++) {
            Pr[i] = 2 * Pc[i] - Pworst[i];
        }
      
        if (f(Pr) < f(Plworst))
        {
           
            for (int i = 0; i < n; i++) {
                Pe[i] = Pc[i] + alpha * (Pr[i] - Pc[i]);

            }
            if (f(Pe) < f(Pr)) {
                for (int i = 0; i < n; i++) {
                    xx[v[n].index] = Pe;
                }
            }
            else {
                for (int i = 0; i < n; i++) {
                    xx[v[n].index] = Pr;
                }
            }
            
        }
        else
        {
           
            for (int i = 0; i < n; i++) {
                Py[i] = Pc[i] + beta * (Pr[i] - Pc[i]);
            }
            
            if (f(Py) < f(Pr))
            {
                xx[v[n].index] = Py;
            }
            else
            {
                
                for (int i = 0; i < n + 1; i++) {
                    for (int j = 0; j < n; j++) {
                        xx[i][j] = Pbest[j] + theta * (xx[i][j] - Pbest[j]);
                    }

                }
             }   
        }
  
    }
}
int main()
{
    int n = 0;
    cin >> n;
    vector<double> x0(n);
    for (int i = 0; i < n; i++) {
        cin >> x0[i];
    }
    double h = 0;
    cin >> h;
    double alpha = 0;
    cin >> alpha;
    double beta = 0;
    cin >> beta;
    double theta = 0;
    cin >> theta;
    double eps = 0;
    cin >> eps;
    vector<double> x_opt(n);
    double f_opt = 0;
    int stepnum = 0;
   
    SimpleSearch(x0, h, n, alpha, beta, theta, eps, x_opt, f_opt, stepnum);
    for (int i = 0; i < n; i++)
    {
        printf("%.6lf ", x_opt[i]);
    }
    printf("%.6lf ", f_opt);
    printf("%d\n", stepnum);
    return 0;
}