用共轭梯度法求函数极小值和最优解，其中用进退法求步长区间，用黄金分割法求最佳步长

//共轭梯度法求二次函数的最优解和极小值/*----------------------------------------------*///本函数为 x1*x1+4*x2*x2-2*x1*x2-4*x1/*----------------------------------------------*///若需要改为其他函数，请修改//change#include<stdio.h>#include<math.h>#include<stdlib.h>#pragma warning(disable:4996)//取消scanf不安全提示//原函数#define f(x1,x2) x1*x1+4*x2*x2-2*x1*x2-4*x1//change//梯度模#define tdm(x1,x2) sqrt((2*x1-2*x2-4)*(2*x1-2*x2-4)+(8*x2-2*x1)*(8*x2-2*x1))//change//x1的偏导数#define G1(x1,x2) 2*x1-2*x2-4//change//x2的偏导数#define G2(x1,x2) 8*x2-2*x1//change//一维搜索//进退法求搜索区间const float eps = 0.00001;//eps为计算精度double HJFC(double x1[], double s1[]){int k = 1, i, j;double a0 = 1, b0 = 0.5, a1, b1, a[3], f[3], y[3][2], m, n, ak2, c;//a0为初始步长，b0为初始步长增量，a1，b1为进退法确定的最终区间；a[0] = a0;a[1] = a0 + b0;for (i = 0; i < 2; i++)for (j = 0; j < 2; j++)y[i][j] = x1[j] + a[i] * s1[j];for (i = 0; i < 2; i++)f[i] = f(y[i][0], y[i][1]);if (f[0]>f[1])while (k){b0 = 2 * b0;a[2] = a[1] + b0;for (j = 0; j < 2; j++)y[2][j] = x1[j] + a[2] * s1[j];f[2] = f(y[2][0], y[2][1]);if (f[2]>f[1]){m = a[0];n = a[2];k = 0;}else{k = 1;a[1] = a[2];f[1] = f[2];}}elsewhile (k){b0 = 2 * b0;a[2] = a[0] - b0;for (j = 0; j < 2; j++)y[2][j] = x1[j] + a[2] * s1[j];f[2] = f(y[2][0], y[2][1]);if (f[2]>f[0]){m = a[2];n = a[1];k = 0;}else{k = 1;a[0] = a[2];f[0] = f[2];}}//黄金分割法求最佳步长a1 = m;b1 = n;a[0] = n - 0.618*(n - m);a[1] = m + 0.618*(n - m);for (i = 0; i < 2; i++)for (j = 0; j < 2; j++)y[i][j] = x1[j] + a[i] * s1[j];for (i = 0; i < 2; i++)f[i] = f(y[i][0], y[i][1]);do{if (f[0] < f[1]){n = a[1];a[1] = a[0];f[1] = f[0];a[0] = n - 0.618*(n - m);for (j = 0; j < 2; j++)y[0][j] = x1[j] + a[0] * s1[j];f[0] = f(y[0][0], y[0][1]);}else{m = a[0];a[0] = a[1];f[0] = f[1];a[1] = m + 0.618*(n - m);for (j = 0; j < 2; j++)y[1][j] = x1[j] + a[1] * s1[j];f[1] = f(y[1][0], y[1][1]);}c = (n - m) / (b1 - a1);} while (fabs(c)>eps);ak2 = (m + n) / 2;return ak2;}//共轭梯度void main(){double x[2], s[2], g[4], bita, arph;//x数组为函数解得转置矩阵；s数组为搜索方向的转置矩阵，g数组为梯度转置矩阵,arph为最优步长;int k = 1;//k为迭代次数;//赋初值printf("请输入初始x1,x2（本函数是：1 1）:\n\n");scanf("%d %d", &x[0], &x[1]);printf("过程如下：\n\n");g[0] = G1(x[0], x[1]);g[1] = G2(x[0], x[1]);while (tdm(x[0], x[1]) > eps)//迭代终止准则{if (k == 1){s[0] = -g[0];s[1] = -g[1];bita = 0;}else{bita = (g[0] * g[0] + g[1] * g[1]) / (g[2] * g[2] + g[3] * g[3]);s[0] = -g[0] + bita*s[0];s[1] = -g[1] + bita*s[1];}arph = HJFC(x, s);x[0] = x[0] + arph*s[0];x[1] = x[1] + arph*s[1];g[2] = g[0];g[3] = g[1];g[0] = G1(x[0], x[1]);g[1] = G2(x[0], x[1]);printf("第%d次迭代:\n\n", k);printf("步长steplength=%f\t", arph);printf("bb=%f\t\n", bita);printf("x[0]=%f\tx[1]=%f\n\n\n", x[0], x[1]);k++;}k--;printf("最后结果:\n");printf("最优解为：\nx[0]=%f\nx[1]=%f\n 最小值:min=%f\n 迭代次数：n=%d\n", x[0], x[1], f(x[0], x[1]), k);}