最速下降法的C语言实现

对于无约束优化问题，一般是采用迭代法进行计算，其标准格式为：

x(k+1)=x(k)+a*s 其中a称作步长，s称作方向。步长a一般可以通过一维不精确线搜索（Armijo准则）计算，而根据方向s选择的不同，无约束优化问题一般有最速下降法、BFGS、共轭梯度法，牛顿及拟牛顿法等等，今天只讲最速下降法。

最速下降法的前进方向是目标函数f(x)的负梯度方向。其C语言代码如下：

// Zhuxu12-4.cpp : 定义控制台应用程序的入口点。

#include "stdafx.h"
#include <math.h>


float fun(float x, float y)
{
    return (1 - x)*(1 - x) + 100 * (y - x*x)*(y - x*x);
}


void gfun(float x, float y, float& g1, float& g2)
{
    g1 = (-2)*(1 - x) - 400*x*(y - x*x);
    g2 = 200*(y - x*x);
}


float min(float a, float b)
{
    return (a<b) ? a : b;
}


int armijo(float x, float y, float rho, float sigma, float dx, float dy)
{
    float x1, y1;
    float g1, g2;
    gfun(x, y, g1, g2);
    int i = 0;
    int imax = 100;
    while (i<=imax)
    {
        x1 = x + pow(rho,i)*dx;
        y1 = y + pow(rho, i)*dy;
        if (fun(x1, y1) <= fun(x, y) + sigma*pow(rho, i)*(g1*dx + g2*dy))
            break;
        i = i + 1;
    }
    return i;
}


int main(int argc, char* argv[])
{
    float x0, y0;
    float x, y;
    printf("Please input the initial value:x0,y0\n");
    scanf_s("%f ,%f", &x0, &y0);
    printf("x=%f,y=%f\n", x0, y0);
    int max_iter = 50000;
    float g1,g2,s1,s2;
    x = x0;
    y = y0;
    float rho=0.5f;
    float sigma = 0.2f;
    float EPS = 1e-4f;
    int i = 0;
    int j = 0;
    while (i <= max_iter)
    {
        i = i + 1;
        printf("第%d次：\n", I);
        gfun(x, y, g1, g2);
        if (sqrt(pow(g1, 2) + pow(g2, 2)) < EPS)
            break;
        s1 = -g1;
        s2 = -g2;
        j= armijo(x, y, rho, sigma, s1, s2);
        printf("搜索次数：%d\t", j);
        printf("最佳步长：%f\n", pow(rho,j));
        x =  x + pow(rho, j)*s1;
        y = y + pow(rho, j)*s2;
    }
    printf("最优解是：x=%f,y=%f\n", x,y);
    return 0;
}