Geatpy是由华南理工大学、华南农业大学、德州奥斯汀公立大学本硕博联合团队开发的一款Python上的遗传和进化算法高性能权威框架。其性能远高于matlab遗传算法工具箱及类似的诸如gatbx、GEATbx、gaot等第三方工具箱。尤其适合需要应用遗传或其他进化算法求解建模问题、与当前的其他项目结合、抑或是需要一个平台来从事进化算法研究与实验的用户。

在前面几篇文章中，我们已经介绍了高性能Python遗传和进化算法框架——Geatpy的使用及一些案例。

Geatpy更新方法，以管理员方式运行：：

pip install --user --upgrade geatpy

若在更新过程中遇到”utf8 decode”错误的问题，是windows下用pip进行安装时遇到的常见问题之一解决方法有很多。可以以管理员方式运行cmd并进行更新。

https://blog.csdn.net/qq_33353186/article/details/82014986

https://blog.csdn.net/qq_33353186/article/details/82020507

https://blog.csdn.net/qq_33353186/article/details/82021750

https://blog.csdn.net/qq_33353186/article/details/82047692

详细的Geatpy官方教程文档在github可以查看：

https://github.com/geatpy-dev/geatpy/tree/master/geatpy/doc/tutorial

也可以到官网上查看：http://www.geatpy.com

================分割线==================

下面就等式约束问题展开论述：

含等式约束的目标优化问题是常见的比较难用遗传算法求解的一类问题。经典的处理方法有以下三类：

一. 直接采取罚函数来处理等式约束：

这种方法相对简单，直接把遗传算法中不满足等式约束条件的个体赋予最低的适应度，使其在种群进化中逐渐被淘汰。

然而，罚函数在处理不等式约束较为成功，但面对等式约束这类可行解空间远比非可行解空间小的问题时，无论罚函数对处于非可行解的个体加以多么严厉的惩罚，种群进化过程中会很难跳出非可行解区域，这就极有可能导致遗传算法出现连可行解都找不到的结果。

二. 将等式约束化成一个新的目标，用多目标优化算法去处理等式约束问题：

例如一个等式约束为 x + y + z = 1，那么，可以设计一个新的目标函数：

该目标表示 x + y + z – 1的绝对值越小越好。该目标越接近0，说明越接近满足约束条件。当该目标的值为0时，对应的解就能满足约束条件了。

这种方法可以极大地拓宽所研究问题的可行解的区域，从而使得遗传算法更方便地进行最优解的搜索。

但是：这种方法也存在弊端。用这种方法对可行解的拓宽程度相对较大，并且，由于使用多目标优化算法求解，它得到的结果是一个帕累托最优解集。而问题最容易出在：有时候会因为在搜索时解的分布性不好而导致的搜索得到的帕累托前沿种并不包含使得上述新增目标达到绝对最优（即abs(x+y+z-1)=0）的解，从而导致所有帕累托最优解均不满足该等式约束。

因此方法二也会导致算法的不稳定。

三. 采用数学方法降维消除等式约束：

该方法是利用数学推导，利用等式约束，将其中一个或多个变量用其余变量线性表出或非线性表出，从而使模型降维，同时也成功消除了等式约束，或是成功地将等式约束化为另一类不等式约束。

这种方法在效果上是最好的。但其局限性是：当对于复杂的等式约束、无法将其中一个或多个变量用其余变量线性表出或非线性表出时，该方法无法使用。

因此，应优先使用方法三，当不可行时再考虑方法一和二。

================分割线==================

采用Geatpy来求解这类含等式约束优化问题：

问题抽象成模型描述如下：

我们利用x3 = 1 – x1 – x2实现降维。此时成功将等式约束条件转化成了x1 + x2 <= 1这个不等式约束条件，此时由于x3 = 1 – x1 – x2 必定小于1，因此降维后模型是等价于原模型的。

Geatpy解决这类约束优化问题极为简单。

（Geatpy v1.1.0版本之后，对数据结构和算法模板作出了重大升级，因此如果使用内置的算法模板，则目标函数aimfuc和适应度罚函数punishing需要按照下面格式定义：）

首先创建目标函数，写在“aimfuc.py”文件中：

# -*- coding: utf-8 -*-
"""
aimfc.py - 目标函数demo
描述:
    Geatpy的目标函数遵循本案例的定义方法，
    若要改变目标函数的输入参数、输出参数的格式，则需要修改或自定义算法模板
"""

def aimfuc(Phen, LegV):
    x1 = Phen[:, [0]]
    x2 = Phen[:, [1]]
    x3 = 1 - x1 - x2 # 将x1 + x2 + x3 = 1的等式约束降维化处理
    f = 4 * x1 + 2 * x2 + x3
    # 约束条件
    idx1 = np.where(2 * x1 + x2 > 1)[0]
    idx2 = np.where(x1 + 2 * x3 > 2)[0]
    idx3 = np.where(x1 + x2 > 1)[0]
    exIdx = np.unique(np.hstack([idx1, idx2, idx3])) # 得到非可行解在种群中的下标
    
    # 惩罚方法1：修改非可行个体的目标函数值的方法，对非可行解对应的目标函数值作出惩罚，而不标记其为非可行解
    # 不标记非可行解的做法比较危险，假如处理不慎，可能会导致以下错误情况：
    # 某一代的种群全是非可行解而其最优的目标函数值恰好是历代最优的情况出现。
    # 此时进化记录器会受到欺骗，而得出一个实际上是非可行解的搜索结果。
    # 一个比较好的处理方法是：要让其绝对地比其他任何可行解的目标函数值都要大或小（看是最小化目标还是最大化目标）。
#    f[exIdx] = -1
    
    # 惩罚方法2： 标记非可行解在可行性列向量中对应的值为0，并编写punishing罚函数来修改非可行解的适应度。
    # 也可以不写punishing，因为Geatpy内置的算法模板及内核已经对LegV标记为0的个体的适应度作出了修改。
    # 使用punishing罚函数实质上是对非可行解个体的适应度作进一步的修改
    LegV[exIdx] = 0 # 对非可行解作出标记，使其在可行性列向量中对应的值为0
    
    return [f,LegV]

编写适应度罚函数punishing：

# -*- coding: utf-8 -*-
"""
punishing.py - 罚函数demo
描述:
    Geatpy的罚函数遵循本案例的定义方法，
    若要改变罚函数的输入参数、输出参数的格式，则需要修改或自定义算法模板
"""

import numpy as np

def punishing(LegV, FitnV):
    FitnV[np.where(LegV == 0)[0]] = 0
    return FitnV

最后编写执行脚本”main.py”：

# -*- coding: utf-8 -*-
"""
执行脚本main.py
描述：
    该demo是展示如何计算带等式约束的单目标优化问题:
        max 4 * x1 + 2 * x2 + x3
        s.t.
            2 * x1 + x2 <= 1
            x1 + 2 * x3 > 2
            x1 + x2 + x3 == 1
            0 <= x1 <= 1
            0 <= x2 <= 1
            0 <= x3 <= 2
    其中目标函数和约束条件写在aimfuc.py文件中，适应度罚函数写在罚函数文件punishing.py中
    本案例通过降维的方法，将等式约束化成了不等式约束，大大拓宽了可行解的空间，方便遗传算法求解
    此外，本案例展示了利用多种群竞争的进化算法模板sga_mpc_real_templet了解决该问题。
"""

import numpy as np
import geatpy as ga

# 获取函数接口地址
AIM_M = __import__('aimfuc')
PUN_M = __import__('punishing')
# 变量设置
x1 = [0, 1] # 自变量1的范围
x2 = [0, 1] # 自变量2的范围
b1 = [1, 1] # 自变量1是否包含下界
b2 = [1, 1] # 自变量2是否包含上界
ranges=np.vstack([x1, x2]).T # 生成自变量的范围矩阵
borders = np.vstack([b1, b2]).T # 生成自变量的边界矩阵
precisions = [2] * 2 # 自变量的编码精度，由于控制变量是连续型的，并且算法模板中使用实值编码，因此这里设置成大于0的任意值即可
newRanges = ga.meshrng(ranges, gridnum = 2) # 对控制变量范围进行网格化，网格边长为2
# 生成网格化后的区域描述器集合
FieldDRs = []
for i in range(len(newRanges)):
    FieldDRs.append(ga.crtfld(newRanges[i], borders, precisions))
# 调用编程模板
[pop_trace, var_trace, times] = ga.sga_mpc_real_templet(AIM_M, 'aimfuc', PUN_M,\
 'punishing', FieldDRs, problem = 'R', maxormin = -1, MAXGEN = 50, NIND = 50,\
 SUBPOP = 1, GGAP = 0.9, selectStyle = 'etour', recombinStyle = 'xovdprs',\
 recopt = 0.9, pm = 0.3, distribute = True, drawing = 1)

执行”main.py”,得到结果如下：

最优的目标函数值为： 2.5
最优的控制变量值为：
0.5
0.0
最优的一代是第 3 代
时间已过 0.4036142826080322 秒

分析：sga_mpc_real_templet源代码见：

https://github.com/geatpy-dev/geatpy/blob/master/geatpy/source-code/templets/sga_mpc_real_templet.py

该模板实现了基于多种群竞争进化的单目标优化算法。由于多种群是在相互竞争中进化的，因此收敛速度极快。并且，每个种群的最大进化代数也不需要设置太大，本例中设置50代就够了。

Geatpy另外还提供了多种群独立进化的单目标优化算法模板:sga_mps_real_templet。与上面的不同的是该算法模板中种群是互相独立进化，最后再挑选最优个体的，而不是在竞争中进化。

在调用多种群进化算法模板sga_mpc_real_templet和sga_mps_real_templet前，常根据各控制变量的范围先对其进行网格化，使其网格化成若干个小的范围，然后不同范围用不同的种群去进化搜索。变量范围网格化的库函数为：meshrng，用法详见：https://github.com/geatpy-dev/geatpy/blob/master/geatpy/doc/API/meshrng/meshrng.pdf

更多的demo详见：

https://github.com/geatpy-dev/geatpy/tree/master/geatpy/demo

欢迎继续跟进，感谢！