//************************************************************************
//**************************共轭梯度算法C程序******************************
//************************************************************************
//***********本程序适用于n设计变量的函数优化问题,对于不同的设计变量*********
//***********个数可以改变维数,该算法程序只在主函数中与其他无约束*************
//***********程序有差别,其他部分基本一样。 ***********
//************************************************************************
//************************************************************************
//************************************************************************
#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#include<conio.h>
#define tt 0.01 //-------------------一维搜索初始步长---------------------------
#define ff 1.0e-6 //-------------------差分法求梯度时的步长-----------------------
#define ac 1.0e-6 //-------------------终止迭代梯度模收敛精度---------------------
#define ad 1.0e-6 //-------------------一维搜索收敛精度---------------------------
#define ax 1.0e-6 //-------------------终止迭代点距收敛精度-----------------------
#define n 2 //-------------------设计变量的维数-----------------------------
double ia;
//=================================输入目标函数表达式========================================
double fny(double *x)
{
double x1=x[0],x2=x[1],x3=x[2];
double f;
f=x1*x1*x1*x1+x1*x2+(1+x2)*(1+x2);
return f;
}
//==========================================================================================
//==================================迭代更新最优点===========================================
double * iterate(double *x,double a,double *s)
{
double *x1;
int i;
x1=(double *)malloc(n*sizeof(double));
for(i=0;i<n;i++)
x1[i]=x[i]+a*s[i];
return x1;
}
//==========================================================================================
//===================================计算更新后节点函数值====================================
double func(double *x,double a,double *s)
{
double *x1;
double f;
x1=iterate(x,a,s);
f=fny(x1);
return f;
}
//==========================================================================================
//===================================确定初始单谷区间的退进法================================
void finding(double a[3],double f[3],double *xk,double *s)
{
double t=tt;
int i;
double a1,f1;
a[0]=0;f[0]=func(xk,a[0],s);
for(i=0;;i++) //循环次数i,受break执行的影响
{
a[1]=a[0]+t;
f[1]=func(xk,a[1],s);
if(f[1]<f[0]) break;
if(fabs(f[1]-f[0])>=ad)
{
t=-t;
a[0]=a[1];f[0]=f[1];
}
else
{
if(ia==1) return; //break
t=t/2;ia=1;
}
}
for(i=0;;i++)
{
a[2]=a[1]+t;
f[2]=func(xk,a[2],s);
if(f[2]>f[1]) break;
t=2*t;
a[0]=a[1];f[0]=f[1];
a[1]=a[2];f[1]=f[2];
}
if(a[0]>a[2])
{
a1=a[0];
f1=f[0];
a[0]=a[2];
f[0]=f[2];
a[2]=a1;
f[2]=f1;
}
return;
}
//==========================================================================================
//=======================================一维搜索============================================
double lagrange(double *xk,double *ft,double *s)
{
int i;
double a[3],f[3];
double b,c,d,aa;
finding(a,f,xk,s);
for(i=0;;i++)
{
if(ia==1) { aa=a[1]; *ft=f[1]; break; }
d=(pow(a[0],2)-pow(a[2],2))*(a[0]-a[1])-(pow(a[0],2)-pow(a[1],2))*(a[0]-a[2]);
if(fabs(d)==0) break;
c=((f[0]-f[2])*(a[0]-a[1])-(f[0]-f[1])*(a[0]-a[2]))/d;
if(fabs(c)==0) break;
b=((f[0]-f[1])-c*(pow(a[0],2)-pow(a[1],2)))/(a[0]-a[1]);
aa=-b/(2*c);
*ft=func(xk,aa,s);
if(fabs(aa-a[1])<=ad) {if(*ft>f[1]) aa=a[1];break;}
if(aa>a[1])
{
if(*ft>f[1]) {a[2]=aa;f[2]=*ft;}
else if(*ft<f[1]) {a[0]=a[1];a[1]=aa;f[0]=f[1];f[1]=*ft;}
else if(*ft==f[1])
{
a[2]=aa;a[0]=a[1];
f[2]=*ft;f[0]=f[1];
a[1]=(a[0]+a[2])/2;
f[1]=func(xk,a[1],s);
}
}
else
{
if(*ft>f[1]) {a[0]=aa;f[0]=*ft;}
else if(*ft<f[1]) {a[2]=a[1];a[1]=aa;f[2]=f[1];f[1]=*ft;}
else if(*ft==f[1])
{a[0]=aa;a[2]=a[1];
f[0]=*ft;f[2]=f[1];
a[1]=(a[0]+a[2])/2;
f[1]=func(xk,a[1],s);
}
}
}
if(*ft>f[1]) {*ft=f[1];aa=a[1];}
return aa;
}
//==========================================================================================
//====================================计算向量梯度===========================================
double *gradient(double *xk)
{
double *g,f1,f2,q;
int i;
g=(double*)malloc(n*sizeof(double));
f1=fny(xk);
for(i=0;i<n;i++)
{q=ff;
xk[i]=xk[i]+q; f2=fny(xk);
g[i]=(f2-f1)/q; xk[i]=xk[i]-q;
}
return g;
}
//==========================================================================================
//=======================================F_R算法主程序=======================================================
double * fr(double *xk)
{
double mod2,mod1,s[n];
double aa,ib;
double *ft,*xk1,*g1,*g2,*xx,*x0=xk;
double fi,fxk1;
int i,j,k,step=1;
ft=(double *)malloc(sizeof(double));
xk1=(double *)malloc(n*sizeof(double));
//=========================================迭代求解最优点=====================================
for(k=0;k<n;k++)
{
if(ia==1) { xx=xk; break; }
if(step==1) //====================求初始点X(0)==============================
{
g1=gradient(xk);
fxk1=fny(xk);
x0=xk;
xk1=xk;
for(i=0;i<n;i++)
s[i]=-g1[i];
} //============================================================
else if(step==2) //=====================求点X(1)=================================
{
aa=lagrange(xk,ft,s);
xk1=iterate(xk,aa,s);
g2=gradient(xk1);
fxk1=fny(xk1);
} //============================================================
else //=====================求X(k+1)点============================
{
mod1=0.0;
mod2=0.0;
for(i=0;i<n;i++)
mod1=mod1+g1[i]*g1[i];
for(i=0;i<n;i++)
mod2=mod2+g2[i]*g2[i];
for(i=0;i<n;i++)
s[i]=-g2[i]+mod2*100.0/mod1/100.0*s[i];
g1=g2;
xk=xk1;
aa=lagrange(xk,ft,s);
xk1=iterate(xk,aa,s);
g2=gradient(xk1);
fxk1=fny(xk1);
} //===========================================================
printf("%d %f %f %f \n",step,xk1[0],xk1[1],fxk1);
if(step>=2)
{
if(sqrt((g2[1]-g1[1])*(g2[1]-g1[1])+(g2[0]-g1[0])*(g2[0]-g1[0]))<ac)
break;
}
step=step+1;
k=-1;
}
//==========================================================================================
xk=x0;
xx=xk1;
return xx;
}
//==========================================================================================================
main () //=========================================主函数=====================================
{
printf("THIS IS THE F_R METHOD !\n\n\n");
int k;
double *xx,f;
double xk[n]={0,0};
xx=fr(xk);
f=fny(xx);
printf("\n\nThe Optimal Design Result Is:\n");
for(k=0;k<n;k++)
{printf("\n\tx[%d]*=%f",k+1,xx[k]);}
printf("\n\tf*=%f",f);
getch();
}
//====================================================================================
//问题:f(X)=x1*x1*x1*x1+x2*x2+1时
//出现程序错误,估计是因为迭代公式中
//出现求梯度模的分量,有待继续分析。
//
gongetidu.rar_共轭梯度_共轭梯度法
版权申诉
124 浏览量
2022-09-20
17:31:55
上传
评论
收藏 2KB RAR 举报
周楷雯
- 粉丝: 80
- 资源: 1万+
最新资源
- JAVA实现Modbus RTU或Modbus TCPIP案例.zip
- 基于YOLOv8的FPS TPS AI自动锁定源码+使用步骤说明.zip
- JAVA实现Modbus RTU或Modbus TCPIP案例.zip
- 基于yolov8+streamlit的火灾检测部署源码+模型.zip
- 测试aaaaaaabbbbb
- VID20240521070643.mp4
- Android系统原理与开发学习要点详解-培训课件.zip
- 部署yolov8的tensorrt模型支持检测分割姿态估计的C++源码+部署步骤.zip
- 以简单、易用、高性能为目标、开源的时序数据库,支持Linux及Windows, Time Series Database.zip
- python-leetcode面试题解之第198题打家劫舍-题解.zip
资源上传下载、课程学习等过程中有任何疑问或建议,欢迎提出宝贵意见哦~我们会及时处理!
点击此处反馈
评论0