遗传算法的C语言实现(一):以非线性函数求极值为例

时间 2019-12-13

标签遗传算法 c语言实现非线性函数极值为例繁體版

原文原文链接

之前搞数学建模的时候，研究过(其实也不算是研究，只是大概了解)一些人工智能算法，好比前面已经说过的粒子群算法(PSO),还有著名的遗传算法(GA),模拟退火算法(SA),蚁群算法(ACA）等。当时懂得很是浅，只会copy别人的代码(通常是MATLAB),改一改值和参数，东拼西凑就拿过来用了，根本没有搞懂的其内在的原理究竟是什么。这一段时间，我又从新翻了一下当时买的那本《MATLAB智能算法30个案例分析》，重读一遍，发现这本书讲的仍是很是不错的，不只有现成算法的MATLAB实现，并且把每一种算法的实现原理都基本讲清楚了，我再次看的时候，之前不懂的(也许之前太急功近利了)如今竟然基本所有弄懂了。真是让我感到很是开心~~~因此我萌发了一个冲动，想用如今本身比较熟悉的C语言，把这本书上的案例从新实现一遍，完全弄懂每一个算法的原理。我打算不按期地更新几篇来介绍其中几种常见的智能算法和其应用(固然我也不止会看这本书，还会参考别的书籍和论文),尽可能把每种算法的原理给讲明白，也算是对我接触这些智能算法的一个总结吧。html

今天首先介绍遗传算法(genetic algorithm,GA)。算法

遗传算法是一种进化算法，其基本原理是模仿天然界中的生物“物竞天择，适者生存”的进化法则，把问题参数编码为染色体，再利用迭代的方式进行选择、交叉、变异等运算法则来交换种群中染色体的信息，最终生成符合优化目标的染色体。ubuntu

在遗传算法中，染色体对应的是数据或者数组，一般是由一维的串结构数据来表示，串上各个位置对应基因的的取值。基因组成的串就是染色体，或者称为基因型个体。必定数量的个体组成了种群，种群中个体数目的大小称为种群大小，也称为种群规模，而各个个体对环境的适应程度称为适应度。数组

　　标准遗传算法的步骤以下：函数

　　(1) 编码：遗传算法在搜索解空间以前须要将解数据表示成遗传空间的基因型串结构数据，这些串结构数据的不一样组合构成了不一样的染色体。优化

　　(2)初始化：即生成初始种群。具体的作法是随机生成N个初始的染色体(解空间的解),每个染色体其实就至关于一个个体，N个个体构成了一个初始种群。遗传算法以这N个个体做为初始值开始进化。编码

　　(3) 适应度评估:适应度代表个体或者解的优劣性。不一样的问题，定义适应度函数的方式也不一样。好比若是是求一个函数的最大值，那么通常某个解对应的函数值越大，那么这个个体(解)的适应度也就越高。人工智能

（4）选择:选择的目的是为了从当前种群中选出优良的个体，使它们有机会成为父代繁殖下一代子孙。遗传算法经过选择过程体现这一思想，进行选择的原则是适应性强的个体为下一代贡献一个或者多个后代的几率大。这体现了达尔文的适者生存原则。spa

　　(5)交叉：交叉操做是遗传算法中最主要的遗传操做。经过交叉能够获得新一代个体，新个体组合了其父代的个体特征。交叉体现了信息交换的思想。code

　　(6)变异:变异首先在群体中随机选择一个个体，对于选中的个体以必定的几率（一般是比较小的几率，这与天然界一致，天然界的变异都是小几率事件）随机改变染色体中某个基因的值。

有的时候除了选择选择、交叉、变异这三种操做以外，咱们还会针对具体的问题加入其它的操做（好比逆转之类），可是选择、交叉、变异是全部的遗传算法都共同的拥有的遗传操做。

本文以遗传算法常见的一个应用领域——求解复杂的非线性函数极值问题为例，来讲明如何用代码(这里是用C语言，固然你能够换作别的任何一种语言)来具体实现上述的遗传算法操做。求解的函数来自《MATLAB智能算法30个案例分析》，即：

f = -5*sin(x1)*sin(x2)*sin(x3)*sin(x4)*sin(x5) - sin(5*x1)*sin(5*x2)*sin(5*x3)*sin(5*x4)*sin(5*x5) + 8

其中x1,x2,x3,x4,x5是0~0.9*PI之间的实数。该函数的最小值为2,当x1,x2,x3,x4,x5都取PI/2时获得。

使用C语言实现的遗传算法求解代码以下：

/*
 * 遗传算法求解函数的最优值问题
 * 参考自《MATLAB智能算法30个案例分析》
 * 本例的寻优函数为:f = -5*sin(x1)*sin(x2)*sin(x3)*sin(x4)*sin(x5) - sin(5*x1)*sin(5*x2)*sin(5*x3)*sin(5*x4)*sin(5*x5) + 8
 * 其中x1,x2,x3,x4,x5是0~0.9*PI之间的实数。该函数的最小值为2,当x1,x2,x3,x4,x5都取PI/2时获得。
 * update in 16/12/3
 * author:Lyrichu
 * email:919987476@qq.com
 */
#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#include<time.h>
#define PI 3.1415926 //圆周率
#define sizepop 50 // 种群数目
#define maxgen 500 // 进化代数
#define pcross 0.6 // 交叉几率
#define pmutation 0.1 // 变异几率
#define lenchrom 5 // 染色体长度
#define bound_down 0 // 变量下界，这里由于都相同，因此就用单个值去代替了，若是每一个变量上下界不一样，也许须要用数组定义
#define bound_up (0.9*3.1415926) // 上界
double chrom[sizepop][lenchrom]; // 种群数组
double fitness[sizepop]; //种群每一个个体的适应度
double fitness_prob[sizepop]; // 每一个个体被选中的几率(使用轮盘赌法肯定)
double bestfitness[maxgen]; // 每一代最优值
double gbest_pos[lenchrom]; // 取最优值的位置
double average_best[maxgen+1]; // 每一代平均最优值
double gbest; // 全部进化中的最优值
int gbest_index; // 取得最优值的迭代次数序号


// 目标函数
double fit_func(double * arr)
{
    double x1 = *arr;
    double x2 = *(arr+1);
    double x3 = *(arr+2);
    double x4 = *(arr+3);
    double x5 = *(arr+4);
    double func_result = -5*sin(x1)*sin(x2)*sin(x3)*sin(x4)*sin(x5) - sin(5*x1)*sin(5*x2)*sin(5*x3)*sin(5*x4)*sin(5*x5) + 8;
    return func_result;
}

// 求和函数
double sum(double * fitness)
{
    double sum_fit = 0;
    for(int i=0;i<sizepop;i++)
        sum_fit += *(fitness+i);
    return sum_fit;
}
// 求最小值函数
double * min(double * fitness)
{
    double min_fit = *fitness;
    double min_index = 0;
    static double arr[2];
    for(int i=1;i<sizepop;i++)
    {
        if(*(fitness+i) < min_fit)
        {
            min_fit = *(fitness+i);
            min_index = i;
        }
    }
    arr[0] = min_index;
    arr[1] = min_fit;
    return arr;
}
// 种群初始化
void init_chrom()
{
    for(int i=0;i<sizepop;i++)
    {
        for(int j=0;j<lenchrom;j++)
        {
            chrom[i][j] = bound_up*(((double)rand())/RAND_MAX);
        }
        fitness[i] = fit_func(chrom[i]); // 初始化适应度
    }
}
// 选择操做
void Select(double chrom[sizepop][lenchrom])
{
    int index[sizepop];
    for(int i=0;i<sizepop;i++)
    {
        double * arr = chrom[i];
        fitness[i] = 1/(fit_func(arr)); // 本例是求最小值，适应度为目标函数的倒数，即函数值越小，适应度越大
    }
    double sum_fitness = 0;
    for(int i=0;i<sizepop;i++)
    {
        sum_fitness += fitness[i]; // 适应度求和
    }
    for(int i=0;i<sizepop;i++)
    {
        fitness_prob[i] = fitness[i]/sum_fitness;
    }
    for(int i=0;i<sizepop;i++)
    {
        fitness[i] = 1/fitness[i]; // 恢复函数值
    }
    for(int i=0;i<sizepop;i++) // sizepop 次轮盘赌
    {
        double pick = ((double)rand())/RAND_MAX;
        while(pick < 0.0001)
            pick = ((double)rand())/RAND_MAX;
        for(int j=0;j<sizepop;j++)
        {
            pick = pick - fitness_prob[j];
            if(pick<=0)
            {
                index[i] = j;
                break;
            }

        }
        
    }
    for(int i=0;i<sizepop;i++)
    {
        for(int j=0;j<lenchrom;j++)
        {
            chrom[i][j] = chrom[index[i]][j];
        }
        fitness[i] = fitness[index[i]];
    }
}

// 交叉操做
void Cross(double chrom[sizepop][lenchrom])
{
    for(int i=0;i<sizepop;i++)
    {
        // 随机选择两个染色体进行交叉
        double pick1 = ((double)rand())/RAND_MAX;
        double pick2 = ((double)rand())/RAND_MAX;
        int choice1 = (int)(pick1*sizepop);// 第一个随机选取的染色体序号
        int choice2 = (int)(pick2*sizepop);// 第二个染色体序号
        while(choice1 > sizepop-1)
        {
            pick1 = ((double)rand())/RAND_MAX;
            choice1 = (int)(pick1*sizepop); //防止选取位置过界
        }
        while(choice2 > sizepop-1)
        {
            pick2 = ((double)rand())/RAND_MAX;
            choice2 = (int)(pick2*sizepop);
        }
        double pick = ((double)rand())/RAND_MAX;// 用于判断是否进行交叉操做
        if(pick>pcross)
            continue;
        int flag = 0; // 判断交叉是否有效的标记
        while(flag == 0)
        {
            double pick = ((double)rand())/RAND_MAX;
            int pos = (int)(pick*lenchrom);
            while(pos > lenchrom-1)
            {
                double pick = ((double)rand())/RAND_MAX;
                int pos = (int)(pick*lenchrom); // 处理越界
            }
            // 交叉开始
            double r = ((double)rand())/RAND_MAX;
            double v1 = chrom[choice1][pos];
            double v2 = chrom[choice2][pos];
            chrom[choice1][pos] = r*v2 + (1-r)*v1;
            chrom[choice2][pos] = r*v1 + (1-r)*v2; // 交叉结束
            if(chrom[choice1][pos] >=bound_down && chrom[choice1][pos]<=bound_up && chrom[choice2][pos] >=bound_down && chrom[choice2][pos]<=bound_up)
                flag = 1; // 交叉有效，退出交叉，不然继续下一次交叉

        }

    }
}

// 变异操做
void Mutation(double chrom[sizepop][lenchrom])
{
    for(int i=0;i<sizepop;i++)
    {
        double pick = ((double)rand())/RAND_MAX; // 随机选择一个染色体进行变异
        int choice = (int)(pick*sizepop);
        while(choice > sizepop-1)
        {
            pick = ((double)rand())/RAND_MAX;
            int choice = (int)(pick*sizepop);  // 处理下标越界
        }
        pick = ((double)rand())/RAND_MAX; // 用于决定该轮是否进行变异
        if(pick>pmutation)
            continue;
        pick = ((double)rand())/RAND_MAX;
        int pos = (int)(pick*lenchrom);
        while(pos > lenchrom-1)
        {
            pick = ((double)rand())/RAND_MAX;
            pos = (int)(pick*lenchrom);
        }
        double v = chrom[i][pos];
        double v1 = v - bound_up;
        double v2 = bound_down - v;
        double r = ((double)rand())/RAND_MAX;
        double r1 = ((double)rand())/RAND_MAX;
        if(r >= 0.5)
            chrom[i][pos] = v - v1*r1*(1-((double)i)/maxgen)*(1-((double)i)/maxgen);
        else
            chrom[i][pos] = v + v2*r1*(1-((double)i)/maxgen)*(1-((double)i)/maxgen); // 注意这里写法是不会越界的，因此只用一次变异就能够了
    }
}

// 主函数
int main(void)
{
    time_t start,finish;
    start = clock(); // 程序开始计时
    srand((unsigned)time(NULL)); // 初始化随机数种子
    init_chrom(); // 初始化种群
    double * best_fit_index = min(fitness);
    int best_index =(int)(*best_fit_index);
    gbest = *(best_fit_index+1); // 最优值
    gbest_index = 0;
    average_best[0] = sum(fitness)/sizepop; //初始平均最优值
    for(int i=0;i<lenchrom;i++)
        gbest_pos[i] = chrom[best_index][i];
    // 进化开始
    for(int i=0;i<maxgen;i++)
    {
        Select(chrom); // 选择
        Cross(chrom); //交叉
        Mutation(chrom); //变异
        for(int j=0;j<sizepop;j++)
        {
            fitness[j] = fit_func(chrom[j]);
        }
        double sum_fit = sum(fitness);
        average_best[i+1] = sum_fit/sizepop; // 平均最优值
        double * arr = min(fitness);
        double new_best = *(arr+1);
        int new_best_index = (int)(*arr); // 新最优值序号
        if(new_best < gbest)
        {
            gbest = new_best;
            for(int j=0;j<lenchrom;j++)
            {
                gbest_pos[j] = chrom[new_best_index][j];
            }
            gbest_index = i+1;
        }

    }
    finish = clock(); // 程序计算结束
    double duration = ((double)(finish-start))/CLOCKS_PER_SEC;
    printf("程序计算耗时:%lf秒\n.",duration);
    printf("遗传算法进化了%d次，最优值为:%lf,最优值在第%d代取得,此代的平均最优值为%lf.\n",maxgen,gbest,gbest_index,average_best[gbest_index]);
    printf("取得最优值的地方为(%lf,%lf,%lf,%lf,%lf).\n",gbest_pos[0],gbest_pos[1],gbest_pos[2],gbest_pos[3],gbest_pos[4]);
    return 0;
}

说明：

（1）这里是求函数的最小值，函数值越小的个体适应度越大，因此这里采用倒数变换，即把适应度函数F(x)定义为1/f(x)(f(x)为目标函数值),这样就保证了适应度函数F(x)越大，个体是越优的。(固然咱们已经保证了f（x）始终是一个正数)

(2) 选择操做采用的是经典的轮盘赌法，即每一个个体被选为父代的几率与其适应度是成正比的，适应度越高，被选中的几率越大，设pi为第i个个体被选中的几率，则

pi=Fi/(∑Fi),Fi为第i个个体的适应度函数。

(3)交叉操做：因为个体采用实数编码，因此交叉操做也采用实数交叉法，第k个染色体a_k和第l个染色体a_l在第j位的交叉操做方法为:

a_kj=a_ij(1-b)+a_ijb; a_lj=a_lj(1-b)+a_kjb,其中b是[0,1]之间的随机数。k，l，j都是随机产生的，即随机选取两个父代，再随机选取一个基因，按照公式完成交叉操做。

(4)变异操做：第i个个体的第j个基因进行变异的操做方法是：

r>=0.5时，a_ij=a_ij+(a_max-a_ij)*f(g)

r<0.5时，a_ij=a_ij+(a_min-a_ij)*f(g)

其中a_max是基因a_ij的上界，a_min是a_ij的下界，f(g)=r2*(1-g/G_max）²，r2是一个[0,1]之间的随机数，g是当前迭代次数，G_max是进化最大次数，r是[0,1]之间的随机数。

能够看出，当r>=0.5时，a_ij是在a_ij到a_max之间变化的；r<0.5时，a_ij是在a_min到a_ij之间变化的。这样的函数其实能够构造不止一种，来完成变异的操做。只要符合随机性，以及保证基因在有效的变化范围内便可。

我在ubuntu16.04下，使用gcc编译器运行获得的结果以下：

从结果上来看，种群数目为50，进化代数为500时，获得的最优解为2.014102已经很是接近实际的最优解2了，固然，若是把种群数目再变大好比500，进化代数变多好比1000代，那么获得的结果将会更精确，种群数目为500，进化代数为1000时的最优解以下：

能够看出，增长种群数目以及增长进化代数以后，最优值从2.014102变为2.000189,比原来的精度提升了两位，可是相应地，程序的运行时间也显著增长了(从0.08秒左右增长到2秒左右)。因此在具体求解复杂问题的时候(好比有的时候会碰到函数参数有几十个之多),咱们就必须综合考虑最优解的精度以及运算复杂度(运算时间)这两个方面，权衡以后，取合适的参数进行问题的求解。

　　固然，这里只是以求解多元非线性函数的极值问题为例来讲明遗传算法的具体实现步骤，遗传算法的实际运用很是普遍，不只仅能够求解复杂函数极值，并且在运筹学等众多领域有着很是普遍的应用，也常常会和其余的智能算法或者优化算法结合来求解问题。这个后面我还会再较为详细地论述。

转自：http://www.cnblogs.com/lyrichu/p/6152897.html