最速下降法 C语言.doc

1.最速下降法#include stdio.h#include math.hdouble fun1(double x1,double x2) /*定义函数fun1为目标函数*/double y; y=x1*x1-2*x1*x2+4*x2*x2+x1-3*x2; return y;void main() double t, x1=1, x2=1,e=0.01, g2, y, m; int k=1; /*定义起始点为x1=0，x2=1，并定义精度为e=0.01*/ g0=2*x1-2*x2+1; /*目标函数对x1求偏导*/ g1=(-2)*x1+8*x2-3; /*目标函数对x2求偏导*/ m=(sqrt(g0*g0+g1*g1); /*对g0*g0+g1*g1求开方，将值赋给m*/ while(me&ke时进行以下循环*/ t=(2*g1-2*g0)*x1+(2*g0-8*g1)*x2-g0+3*g1)/(4*g0*g1-2*g0*g0-8*g1*g1); /*根据梯度法（最速下降法），利用梯度和海赛矩阵*/ x1=x1-g0*t; /*求步长t。 根据步长和梯度方向求出新的x1，x2*/ x2=x2-g1*t;printf(迭代次数%dn,k);printf(搜索方向-%f,-%f,负梯度的模%f,步长%fn,g0,g1,m,t);printf(x的值%f,%fn,x1,x2); g0=2*x1-2*x2+1; g1=(-2)*x1+8*x2-3; m=(sqrt(g0*g0+g1*g1); /*计算新的m*/ printf(新的负梯度的模%fn,m); k+; y=fun1(x1,x2); /*当m不满足me的时候退出循环，并计算fun1,*/ printf(分别输出x1,x2 %f,%fn,x1,x2); /*将值赋给y,并输出。*/ printf(极小值y%f,y);2.FR共轭梯度法#include stdio.h#include math.hdouble fun1(double x1,double x2) /*定义函数fun1为目标函数*/double y; y=x1*x1-2*x1*x2+4*x2*x2+x1-3*x2; return y;double fun2(double g,double d) /*定义函数fun2为求步长的函数*/ double buchang; buchang=-(g0*d0+g1*d1)/(d0*(2*d0-2*d1)+d1*(-2)*d0+8*d1); return buchang;void main() double t, beta, x1=1, x2=1, d2,g4, y, m,e=0.01; /*定义起始点为x1=0，x2=1，并定义精度为e=0.01*/int k=1; g0=2*x1-2*x2+1; /*目标函数对x1求偏导*/ g1=(-2)*x1+8*x2-3; /*目标函数对x2求偏导,求梯度*/ m=(sqrt(g0*g0+g1*g1); /*对g0*g0+g1*g1求开方，将值赋给m*/ while(me&ke时进行以下循环*/ if (k=1) d0=-g0; d1=-g1; beta=0; /*计算因子beta*/elsebeta=(g0*g0+g1*g1)/(g2*g2+g3*g3); /*计算因子beta*/ d0=-g0+beta*d0; d1=-g1+beta*d1; t=fun2(g,d); /*计算步长*/ x1=x1+d0*t; /*根据步长和搜索方向求出新的x1，x2*/ x2=x2+d1*t;printf(迭代次数为%dn,k); printf(梯度%f,%f,梯度的模%fn,g0,g1,m);printf(搜索方向%f,%f,计算因子 %f,步长%fn,d0,d1,beta,t);printf(x的值%f,%fn,x1,x2); g2=g0;g3=g1; g0=2*x1-2*x2+1; /*根据得到的x1，x2求出新的梯度，并将值*/ g1=(-2)*x1+8*x2-3; /*赋给g0，g1,*/ m=double(sqrt(g0*g0+g1*g1); /*计算新的m*/ printf(新的负梯度的模%fn,m);k+; y=fun1(x1,x2); /*当m不满足me的时候退出循环，并计算fun1,*/ printf(分别输出x1,x2 %f,%fn,x1,x2); /*将值赋给y,并输出。*/ printf(极小值y %f,y); 3.DFP法#include#include#include# define e 0.01double fun1 (double x) /*定义目标函数为fun1函数*/ return (double)pow(x0,2)+4*pow(x1,2)-2*x0*x1+x0-3*x1;void fun2 (double x,double g_x) /*定义求梯度函数为fun2函数*/ g_x0 = (double)2*x0-2*x1+1; g_x1 = (double)8*x1-2*x0-3;void fun3 (double h, double g, double d) /*定义求搜索方向函数为fun3函数*/ d0 = -(h0*g0+h1*g1); d1 = -(h2*g0+h3*g1);double fun4 ( double x, double d) /*定义求步长函数为fun4函数*/ return -(2*x0*d0+8*x1*d1-2*x0*d1-2*x1*d0+d0-3*d1) / (2*d0*d0+8*d1*d1-4*d0*d1); void fun5(double h,double x,double tx,double g_x,double tg ) /*定义求H矩阵函数为fun5函数*/ double q2,p2,temp2,th4; q0 = tg0 - g_x0; /*第k+1个点梯度与第k个点处的梯度之差*/ q1 = tg1 - g_x1; p0 = tx0 - x0; /*p(k),第k+1个点与第k个点的点差 */ p1 = tx1 - x1; th0 = h0; /*用临时变量储存原H矩阵的值*/ th1 = h1; th2 = h2; th3 = h3; temp0 = p0 * q0 + p1 * q1; temp1 = (q0 * h0 + q1 * h2) * q0 + (q0 * h1 + q1 * h3) * q1; h0 = th0 + pow(p0,2) / temp0 - pow(th0 * q0 + th1 * q1),2) / temp1; h1 = th1 + (p0 * p1) / temp0 - (th0 * q0 + th1 * q1) * (th2 * q0+ th3 * q1) / temp1; h2 = th2 + (p0 * p1) / temp0 - (th0 * q0 + th1 * q1) * (th2 * q0+ th3 * q1) / temp1; h3 = th3 + pow(p1,2) / temp0 - pow(th2 * q0 + th3 * q1),2) / temp1;void main()int k=1; double x2 = 1 , 1, /*计算初始点*/tx2 = 1 , 1, /*临时存储计算点*/g_x2, /*梯度*/tg2, /*临时梯度*/t, /*步长*/d2, /*搜索方向*/h4 = 1, 0, 0, 1; /*初始H矩阵*/double fx = 0; fun2(x,g_x); /*求梯度*/ tg0=g_x0;tg1=g_x1;fun3(h, g_x, d); /*求搜索方向*/int counter=1; while(sqrt(pow(g_x0,2)+pow(g_x1,2)-e0) do fun3(h, tg, d); /*求搜索方向*/ t = fun4( x, d); /*求步长*/printf(n迭代次数为%dn,counter+);printf(梯度%f,%f,梯度的模%fn,g_x0,g_x1, sqrt(pow(g_x0,2)+pow(g_x1,2); printf(搜索方向%f,%f,步长%fn,d0,d1,t);printf(H矩阵n%f,%fn%f,%f,h0,h1,h2,h3); if(k=2) tx0=x0+t*d0;tx1=x1+t*d1; else tx0=tx0+t*d0; tx1=tx1+t*d1; fun2( tx, tg); /*求梯度*/ fun5(h, x, tx, g_x, tg);/*H矩阵*/if(k = 2) x0 = tx0;x1 = tx1; printf(nx1=%ftx2=%f,x0,x1);fun2( x, g_x);printf(n新梯度的模：%f,sqrt(pow(g_x0,2)+pow(g_x1,2); k+;while(k=2);k = 1; fx = fun1( x);printf(n最优解为：n);printf(nx1=%f, x2=%fn,x0,x1); printf(n最优值为：n); printf(nf(x)=%f n,fx);4.BFGS法#include#include#include# define e 0.01double fun1 (double x) /*定义目标函数为fun1函数*/ return (double)pow(x0,2)+4*pow(x1,2)-2*x0*x1+x0-3*x1;void fun2 (double x,double g_x) /*定义求梯度函数为fun2函数*/ g_x0 = (double)2*x0-2*x1+1; g_x1 = (double)8*x1-2*x0-3;void fun3 (double h, double g, double d) /*定义求搜索方向函数为fun3函数*/ d0 = -(h0*g0+h1*g1); d1 = -(h2*g0+h3*g1);double fun4 ( double x, double d) /*定义求步长函数为fun4函数*/ return -(2*x0*d0+8*x1*d1-2*x0*d1-2*x1*d0+d0-3*d1) / (2*d0*d0+8*d1*d1-4*d0*d1); void fun5(double h,double x,double tx,double g_x,double tg ) /*定义求H矩阵函数为fun5函数*/ double q2,p2,temp2,th4; q0 = tg0 - g_x0; /*第k+1个点梯度与第k个点处的梯度之差*/ q1 = tg1 - g_x1; p0 = tx0 - x0; /*p(k),第k+1个点与第k个点的点差 */ p1 = tx1 - x1; th0 = h0; /*用临时变量储存原H矩阵的值*/ th1 = h1; th2 = h2; th3 = h3; temp0 = p0 * q0 + p1 * q1;/pt*q temp1 = (q0 * h0 + q1 * h2) * q0 + (q0 * h1 + q1 * h3) * q1;/qt*H*q h0 =th0+ (1 + temp1/temp0)*(pow(p0,2) / temp0) - (p0*q0*h0+p0*q1*h2)+(h0*q0*p0+h1*q1*p0) /temp0; h1 =th1+ (1 + temp1/temp0)*(p0 * p1) / temp0) - (p0*q0*h1+p0*q1*h3)+(h0*q0*p1+h1*q1*p1) /temp0; h2 =th2+ (1 + temp1/temp0)*(p1 * p0) / temp0) - (p1*q0*h0+p1*q1*h2)+(h2*q0*p0+h3*q1*p0) /temp0; h3 =th3+ (1 + temp1/temp0)*(p1 * p1) / temp0) - (p1*q0*h1+p1*q1*h3)+(h2*q0*p1+h3*q1*p1) /temp0; void main()int k=1; double x2 = 1 , 1, /*计算初始点*/tx2 = 1 , 1, /*临时存储计算点*/g_x2, /*梯度*/tg2, /*临时梯度*/t, /*步长*/d2, /*搜索方向*/h4 = 1, 0, 0, 1; /*初始H矩阵*/double fx = 0; fun2( x, g_x); /*求梯度*/ tg0=g_x0;tg1=g_x1;fun3(h, g_x, d); /*求搜索方向*/int counter=1;int j=0; while(sqrt(pow(g_x