最优化课程设计
姓 名：指导老师：学 号：班 级：
题目：共轭梯度法
田鑫 智红英 [1**********]6
信息与计算科学111802班
共轭梯度法（Conjugate Gradient）
是介于最速下降法与牛顿法之间的一个方法，它仅需利用一阶导数信息，但克服了最速下降法收敛慢的缺点，又避免了牛顿法需要存储和计算Hesse矩阵并求逆的缺点，共轭梯度法不仅是解决大型线性方程组最有用的方法之一，也是解大型非线性最优化最有效的算法之一。 在各种优化算法中，共轭梯度法是非常重要的一种。其优点是所需存储量小，具有步收敛性，稳定性高，而且不需要任何外来参数。
设我们要求解下列线性系统
其中n-×-n矩阵A是对称的（也即，AT = A），正定的（也即，xTAx > 0对于所有非0向量x属于Rn），并且是实系数的。 将系统的唯一解记作x*。 最后算法
经过一些简化，可以得到下列求解Ax = b的算法，其中A是实对称正定矩阵。
x0 := 0 k := 0 r0 := b
repeat until rk is "sufficiently small":
k := k + 1 if k = 1
p1 := r0 else
end if
xk := xk-1 + αk pk rk := rk-1 - αk A pk end repeat 结果为xk
共轭梯度法程序源代码
#include #include #define N 10
#define eps pow(10,-6)
double f(double x[],double p[],double t) {
double s;
s=pow(x[0]+t*p[0],2)+25*pow(x[1]+t*p[1],2); return s; }
/*以下是进退法搜索区间源程序*/
void sb(double *a,double *b,double x[],double p[]) {
double t0,t1,t,h,alpha,f0,f1; int k=0;
t0=2.5; /*初始值*/ h=1; /*初始步长*/ alpha=2; /*加步系数*/
f0=f(x,p,t0);
f1=f(x,p,t1); while(1) {
if(f1
h=alpha*h; t=t0; t0=t1; f0=f1; k++; } else {
if(k==0) {h=-h;t=t1;} else {
*a=tt1?t:t1; break; } }
t1=t0+h; f1=f(x,p,t1); }
}/*以下是黄金分割法程序源代码*/ double hjfg(double x[],double p[]) {
double beta,t1,t2,t; double f1,f2; double a=0,b=0; double *c,*d; c=&a,d=&b;
sb(c,d,x,p);/*调用进退法搜索区间*/
printf("\nx1=%f,x2=%f,p1=%f,p2=%f",x[0],x[1],p[0],p[1]); printf("\n[a,b]=[%f,%f]",a,b); beta=(sqrt(5)-1.0)/2;
t2=a+beta*(b-a); f2=f(x,p,t2); t1=a+b-t2; f1=f(x,p,t1); while(1) {
if(fabs(t1-t2)
{
t=(t1+t2)/2; b=t2; t2=t1;
f2=f1; t1=a+b-t2; f1=f(x,p,t1); } else {
a=t1; t1=t2; f1=f2;
t2=a+beta*(b-a); f2=f(x,p,t2); } } }
t=(t1+t2)/2; return t; }
/*以下是共轭梯度法程序源代码*/
void gtd() {
double x[N],g[N],p[N],t=0,f0,mod1=0,mod2=0,nanda=0; int i,k,n;
printf("请输入函数的元数值n=2"); scanf("%d",&n);
printf("\n请输入初始值:2,2");
for(i=0;i
g[0]=2*x[0]; g[1]=50*x[1];
mod1=sqrt(pow(g[0],2)+pow(g[1],2));/*求梯度的长度*/ if(mod1>0) {
p[0]=-g[0]; p[1]=-g[1]; k=0; while(1) {
t=hjfg(x,p);/*调用黄金分割法求t的值*/
printf("\np1=%f,p2=%f,t=%f",p[0],p[1],t);
x[0]=x[0]+t*p[0]; x[1]=x[1]+t*p[1]; g[0]=2*x[0]; g[1]=50*x[1];
/*printf("\nx1=%f,x2=%f,g1=%f,g2=%f",x[0],x[1],g[0],g[1]);*/ mod2=sqrt(pow(g[0],2)+pow(g[1],2)); /*求梯度的长度*/
if(mod2
if(k+1==n) {
g[0]=2*x[0]; g[1]=50*x[1]; p[0]=-g[0]; p[1]=-g[1]; k=0; } else {
nanda=pow(mod2,2)/pow(mod1,2);
printf("\nnanda=%f,mod=%f",nanda,mod2); p[0]=-g[0]+nanda*p[0]; p[1]=-g[1]+nanda*p[1]; mod1=mod2; k++; } }
printf("\n--------------------------"); } }
printf("\n最优解为x1=%f,x2=%f",x[0],x[1]); printf("\n最终的函数值为%f",f(x,g,t)); }
main() { gtd(); }
运行结果如下：
请输入函数的元数值n=2 请输入初始值:2 2
x1=2.000000,x2=2.000000,p1=-4.000000,p2=-100.000000 [a,b]=[-4.500000,1.500000]
p1=-4.000000,p2=-100.000000,t=0.020030 nanda=0.001474,mod=3.842730 --------------------------
x1=1.919879,x2=-0.003022,p1=-3.845655,p2=0.003665 [a,b]=[-4.500000,1.500000]
p1=-3.845655,p2=0.003665,t=0.499240 --------------------------
x1=-0.000026,x2=-0.001192,p1=0.000052,p2=0.059610 [a,b]=[-4.500000,1.500000]
p1=0.000052,p2=0.059610,t=0.020000 nanda=0.000000,mod=0.000050 --------------------------
x1=-0.000025,x2=-0.000000,p1=0.000050,p2=0.000001 [a,b]=[-4.500000,1.500000]
p1=0.000050,p2=0.000001,t=0.495505 --------------------------
x1=-0.000000,x2=0.000000,p1=0.000000,p2=-0.000023 [a,b]=[-4.500000,1.500000]
p1=0.000000,p2=-0.000023,t=0.020007 最优解为x1=-0.000000,x2=-0.000000 最终的函数值为0.000000
最优化课程设计
姓 名：指导老师：学 号：班 级：
题目：共轭梯度法
田鑫 智红英 [1**********]6
信息与计算科学111802班
共轭梯度法（Conjugate Gradient）
是介于最速下降法与牛顿法之间的一个方法，它仅需利用一阶导数信息，但克服了最速下降法收敛慢的缺点，又避免了牛顿法需要存储和计算Hesse矩阵并求逆的缺点，共轭梯度法不仅是解决大型线性方程组最有用的方法之一，也是解大型非线性最优化最有效的算法之一。 在各种优化算法中，共轭梯度法是非常重要的一种。其优点是所需存储量小，具有步收敛性，稳定性高，而且不需要任何外来参数。
设我们要求解下列线性系统
其中n-×-n矩阵A是对称的（也即，AT = A），正定的（也即，xTAx > 0对于所有非0向量x属于Rn），并且是实系数的。 将系统的唯一解记作x*。 最后算法
经过一些简化，可以得到下列求解Ax = b的算法，其中A是实对称正定矩阵。
x0 := 0 k := 0 r0 := b
repeat until rk is "sufficiently small":
k := k + 1 if k = 1
p1 := r0 else
end if
xk := xk-1 + αk pk rk := rk-1 - αk A pk end repeat 结果为xk
共轭梯度法程序源代码
#include #include #define N 10
#define eps pow(10,-6)
double f(double x[],double p[],double t) {
double s;
s=pow(x[0]+t*p[0],2)+25*pow(x[1]+t*p[1],2); return s; }
/*以下是进退法搜索区间源程序*/
void sb(double *a,double *b,double x[],double p[]) {
double t0,t1,t,h,alpha,f0,f1; int k=0;
t0=2.5; /*初始值*/ h=1; /*初始步长*/ alpha=2; /*加步系数*/
f0=f(x,p,t0);
f1=f(x,p,t1); while(1) {
if(f1
h=alpha*h; t=t0; t0=t1; f0=f1; k++; } else {
if(k==0) {h=-h;t=t1;} else {
*a=tt1?t:t1; break; } }
t1=t0+h; f1=f(x,p,t1); }
}/*以下是黄金分割法程序源代码*/ double hjfg(double x[],double p[]) {
double beta,t1,t2,t; double f1,f2; double a=0,b=0; double *c,*d; c=&a,d=&b;
sb(c,d,x,p);/*调用进退法搜索区间*/
printf("\nx1=%f,x2=%f,p1=%f,p2=%f",x[0],x[1],p[0],p[1]); printf("\n[a,b]=[%f,%f]",a,b); beta=(sqrt(5)-1.0)/2;
t2=a+beta*(b-a); f2=f(x,p,t2); t1=a+b-t2; f1=f(x,p,t1); while(1) {
if(fabs(t1-t2)
{
t=(t1+t2)/2; b=t2; t2=t1;
f2=f1; t1=a+b-t2; f1=f(x,p,t1); } else {
a=t1; t1=t2; f1=f2;
t2=a+beta*(b-a); f2=f(x,p,t2); } } }
t=(t1+t2)/2; return t; }
/*以下是共轭梯度法程序源代码*/
void gtd() {
double x[N],g[N],p[N],t=0,f0,mod1=0,mod2=0,nanda=0; int i,k,n;
printf("请输入函数的元数值n=2"); scanf("%d",&n);
printf("\n请输入初始值:2,2");
for(i=0;i
g[0]=2*x[0]; g[1]=50*x[1];
mod1=sqrt(pow(g[0],2)+pow(g[1],2));/*求梯度的长度*/ if(mod1>0) {
p[0]=-g[0]; p[1]=-g[1]; k=0; while(1) {
t=hjfg(x,p);/*调用黄金分割法求t的值*/
printf("\np1=%f,p2=%f,t=%f",p[0],p[1],t);
x[0]=x[0]+t*p[0]; x[1]=x[1]+t*p[1]; g[0]=2*x[0]; g[1]=50*x[1];
/*printf("\nx1=%f,x2=%f,g1=%f,g2=%f",x[0],x[1],g[0],g[1]);*/ mod2=sqrt(pow(g[0],2)+pow(g[1],2)); /*求梯度的长度*/
if(mod2
if(k+1==n) {
g[0]=2*x[0]; g[1]=50*x[1]; p[0]=-g[0]; p[1]=-g[1]; k=0; } else {
nanda=pow(mod2,2)/pow(mod1,2);
printf("\nnanda=%f,mod=%f",nanda,mod2); p[0]=-g[0]+nanda*p[0]; p[1]=-g[1]+nanda*p[1]; mod1=mod2; k++; } }
printf("\n--------------------------"); } }
printf("\n最优解为x1=%f,x2=%f",x[0],x[1]); printf("\n最终的函数值为%f",f(x,g,t)); }
main() { gtd(); }
运行结果如下：
请输入函数的元数值n=2 请输入初始值:2 2
x1=2.000000,x2=2.000000,p1=-4.000000,p2=-100.000000 [a,b]=[-4.500000,1.500000]
p1=-4.000000,p2=-100.000000,t=0.020030 nanda=0.001474,mod=3.842730 --------------------------
x1=1.919879,x2=-0.003022,p1=-3.845655,p2=0.003665 [a,b]=[-4.500000,1.500000]
p1=-3.845655,p2=0.003665,t=0.499240 --------------------------
x1=-0.000026,x2=-0.001192,p1=0.000052,p2=0.059610 [a,b]=[-4.500000,1.500000]
p1=0.000052,p2=0.059610,t=0.020000 nanda=0.000000,mod=0.000050 --------------------------
x1=-0.000025,x2=-0.000000,p1=0.000050,p2=0.000001 [a,b]=[-4.500000,1.500000]
p1=0.000050,p2=0.000001,t=0.495505 --------------------------
x1=-0.000000,x2=0.000000,p1=0.000000,p2=-0.000023 [a,b]=[-4.500000,1.500000]
p1=0.000000,p2=-0.000023,t=0.020007 最优解为x1=-0.000000,x2=-0.000000 最终的函数值为0.000000