请选择 进入手机版 | 继续访问电脑版

Hello Mat

 找回密码
 立即注册

QQ登录

只需一步,快速开始

查看: 2729|回复: 1

基于粒子群算法PSO的极小值求解

[复制链接]

806

主题

951

帖子

3万

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
30251
发表于 2017-3-16 23:07:10 | 显示全部楼层 |阅读模式
  1. #include<stdio.h>
  2. #include<string.h>
  3. #include<stdlib.h>
  4. #include<time.h>
  5. #include<math.h>

  6. /* PSO粒子群算法中的两个参数 */
  7. #define c1 2         /* 学习因子1 */
  8. #define c2 2         /* 学习因子2 */
  9. #define w 4          /* 惯性权重  */
  10. #define iter 150     /* 进化次数  */  
  11. #define sizepop 300    /* 种群规模  */
  12. #define Vmax 1       /* 粒子速度上限 */
  13. #define Vmin -1      /* 粒子速度下限 */
  14. #define popmax 9   /* 粒子取值范围上限 */
  15. #define popmin -9  /* 粒子取值范围下限 */

  16. #define length_data 1000   /* 预定义原始待分析数据长度 */

  17. /* 适应度函数 */
  18. double fun(double data[][4],int lineCount,double pop[])  // 染色体的适应度
  19. {
  20.         // 输入:
  21.                // data为原始数据 TOA TR Peo EUI
  22.                // lineCount为数据的长度
  23.                // pop 种群粒子个体
  24.         // 输出:
  25.                // 返回粒子对应的适应度值
  26.         int i;     // 数组长度
  27.         double C0,C1,C2,C3,C4,C5,C6,C7,C8,C9;                   /* 待优化求解未知量 */
  28.         double EUI_pso,TOA,TR,Peo,EUI,Z[length_data],sum_Z=0;   // 原始数据,4组数据
  29.         C0 = pop[0]; C1 = pop[1];  C2 = pop[2]; C3 = pop[3]; C4 = pop[4];
  30.         C5 = pop[5]; C6 = pop[6];  C7 = pop[7]; C8 = pop[8]; C9 = pop[9];
  31. //        leng = sizeof(data)/sizeof(data[0]);                    // 数组的长度--行
  32. //        printf("%lf\t%lf\n",C8,C9);         // 输出显示
  33.         for(i=0;i<lineCount;i++)
  34.         {
  35.                 TOA = data[i][0]; TR = data[i][1]; Peo = data[i][2]; EUI = data[i][3];
  36.         //        printf("%lf\n",TOA);
  37.                 EUI_pso = C0 +C1*TOA +C2*TR +C3*Peo +C4*TOA*TOA +C5*TR*TR +C6*Peo*Peo +C7*TR*Peo +C8*Peo*TOA +C9*TR*TOA;  // PSO算法求解的EUI
  38.                 //printf("%lf\n",EUI_pso);
  39.                 Z[i]= abs( EUI-EUI_pso );  // 绝对误差和
  40.                 // printf("%lf\n",Z[i]);
  41.                 sum_Z = sum_Z + Z[i];      // 叠加和
  42.         }
  43. //        printf("%lf\t\n",sum_Z);         // 输出显示
  44.         
  45.         return sum_Z;                  // 绝对误差和---适应度
  46.         
  47. }

  48. void main(void)
  49. {
  50.         char path_excel[]="F:\\intercompilation\\C\\visual_stUdio_2010\\ysw_vs1\\ysw20140921_c++_3\\PSO3\\PSO3"; //
  51.         FILE *fp;                                     // 文件
  52.         double data[length_data][4],TOA,TR,Peo,EUI;   // 原始数据,4组数据
  53.         int i,j,k,lineCount=0;                        // 上下标
  54.         double rand_a,rand_a1,rand_a2,rand_a3;                   // 随机数
  55.         double pop[sizepop][10],V[sizepop][10],fitness[sizepop]; // 初始化种群变量
  56.         double bestfitness,zbest[10],gbest[sizepop][10],fitnessgbest[sizepop],fitnesszbest,Z;

  57.         // 读取文本文件
  58.         fp = fopen(strcat(path_excel,"\\data.txt"),"r");
  59.         if(fp==NULL)
  60.         {
  61.                 printf("fopen error!\n");
  62.                 return;
  63.         }
  64.         while(!feof(fp))// 判断获取的字符是否是文件的结束符
  65.         {
  66.                 fscanf(fp,"%lf%lf%lf%lf",&TOA,&TR,&Peo,&EUI);
  67.         //        printf("%lf\t%lf\t%lf\t%lf\n",TOA,TR,Peo,EUI);  // 输出显示
  68.                 data[lineCount][0] = TOA;   data[lineCount][1] = TR;   
  69.                 data[lineCount][2] = Peo;   data[lineCount][3] = EUI;
  70.         //        printf("%lf\t%lf\t%lf\t%lf\n",data[lineCount][0],data[lineCount][1],data[lineCount][2],data[lineCount][3]);  // 输出显示
  71.                 lineCount++;
  72.          }
  73.         srand(time(0)); // 初始化随机种子
  74.         // 产生初始粒子和速度
  75.         for(i=0;i<sizepop;i++)
  76.         {
  77.                 // 随机产生一个种群
  78.                 for(j=0;j<10;j++)
  79.                 {
  80.                         rand_a = 2.0*rand()/RAND_MAX-1; // 产生-1到1之间随机数
  81.                 //        printf("%lf\n",rand_a);         // 输出显示
  82.                         pop[i][j] = popmax*rand_a;      // 初始种群   C0-C9十个未知量
  83.                 //        printf("%lf\n",pop[i][j]);
  84.                         V[i][j] = Vmax*rand_a;          // 初始化速度
  85.                 }
  86.         }
  87.         for(i=0;i<sizepop;i++)
  88.         {
  89.                 // 计算适应度
  90.                 fitness[i] = fun(data,lineCount-1,pop[i]);   // 染色体的适应度
  91.         //        printf("%lf\n",fitness[i]);         // 输出显示
  92.         }

  93.         // 找最好的粒子个体--初始化
  94.         bestfitness = fitness[0];         // 绝对误差和最小--初始化
  95.         for(i=0;i<sizepop;i++)
  96.         {
  97.                 for(j=0;j<10;j++)
  98.                 {
  99.                         if(i==0)
  100.                         {
  101.                                 zbest[j] = pop[i][j];   // 全局最佳--初始化
  102.                         }
  103.                         gbest[i][j] = pop[i][j];    // 个体最佳--初始化
  104.                 }
  105.                 fitnessgbest[i] = fitness[i];   // 个体最佳适应度值--初始化
  106.         }
  107.         fitnesszbest = bestfitness;        // 全局最佳适应度值--初始化
  108.         
  109.         // 迭代寻优
  110.         for(i=0;i<iter;i++)
  111.         {
  112.                 for(j=0;j<sizepop;j++)
  113.                 {
  114.                         // 速度更新
  115.                         for(k=0;k<10;k++)
  116.                         {
  117.                                 rand_a1 = rand()/RAND_MAX;   // 产生0到1之间随机数
  118.                                 rand_a2 = rand()/RAND_MAX;   // 产生0到1之间随机数
  119.                                 V[j][k] = w*V[j][k] + c1*rand_a1*(gbest[j][k] - pop[j][k]) + c2*rand_a2*(zbest[k] - pop[j][k]);
  120.                                 if(V[j][k]>Vmax){V[j][k]=Vmax;}  // 速度上限
  121.                                 if(V[j][k]<Vmin){V[j][k]=Vmin;}  // 速度下限
  122.                                 // 种群更新
  123.                                 pop[j][k] = pop[j][k]+0.5*V[j][k];
  124.                                 if(pop[j][k]>popmax){pop[j][k]=popmax;}  // 种群取值上限
  125.                                 if(pop[j][k]<popmin){pop[j][k]=popmin;}  // 种群取值下限
  126.                         }
  127.                         // 自适应变异
  128.                         rand_a3 = rand()/RAND_MAX;   // 产生0到1之间随机数
  129.                         if(rand_a3>0.8)
  130.                         {
  131.                                 k=rand()%10;  // 变异个体标号 0-9
  132.                                 pop[j][k] = popmax*rand()/RAND_MAX;
  133.                         }
  134.                         // 适应度值data,lineCount-1,pop[i])
  135.                         fitness[j]=fun(data,lineCount-1,pop[j]);
  136.                         // 个体最优更新
  137.                         if(fitness[j] < fitnessgbest[j])
  138.                         {
  139.                                 for(k=0;k<10;k++){ gbest[j][k] = pop[j][k]; }
  140.                                 fitnessgbest[j] = fitness[j];
  141.                         }
  142.                         // 群体最优更新
  143.                          if(fitness[j] < fitnesszbest)
  144.                          {
  145.                                  for(k=0;k<10;k++){ zbest[k] = pop[j][k]; }
  146.                                  fitnesszbest = fitness[j];
  147.                          }
  148.                 }  
  149.                 Z = fitnesszbest;  // 存取最好适应度值

  150.         }

  151.         printf("最小绝对误差和值:%lf\n",Z);         // 最好适应度值输出显示
  152.         printf("C0=%lf\tC1=%lf\tC2=%lf\tC3=%lf\tC4=%lf\n",zbest[0],zbest[1],zbest[2],zbest[3],zbest[4]);   // C0-C4输出显示
  153.         printf("C5=%lf\tC6=%lf\tC7=%lf\tC8=%lf\tC9=%lf\n",zbest[5],zbest[6],zbest[7],zbest[8],zbest[9]);   // C5-C9输出显示

  154.         printf("\n");  // 输出显示
  155.         fclose(fp);
  156.         system("pause");    // 消除屏幕一闪即消失的情况
  157. }
复制代码
视频学习:链接:https://pan.baidu.com/s/132jBgG0JLVRAt-b4MM0I9w 密码:pa58
算法QQ  3283892722
群智能算法链接http://halcom.cn/forum.php?mod=forumdisplay&fid=73
回复

使用道具 举报

0

主题

0

帖子

11

积分

新手上路

Rank: 1

积分
11
发表于 2020-5-8 11:38:36 | 显示全部楼层
我的多目标优化结果优于目标:'(
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

Python|Opencv|MATLAB|Halcom.cn  

GMT+8, 2020-6-2 12:16 , Processed in 0.137326 second(s), 23 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.

快速回复 返回顶部 返回列表