粒子群算法尝试解决整数规划问题（效果不好的尝试）

1问题描述

某厂每日8小时的产量不低于1800件．为了进行质量控制，计划聘请两种不同水平的检验员．且每种检验员的日产量不高于1800件．一级检验员的标准为：速度25件/小时，正确率98%，计时工资4元/小时；二级检验员的标准为：速度15件/小时，正确率95%，计时工资3元/小时．检验员每错检一次，工厂要损失2元．为使总检验费用最省，该工厂应聘一级、二级检验员各几名？

2 分析问题

     设需要一级和二级检验员的人数分别是x1 x2，则一天要付给检验员的费用是：

     8*4*x1+8*3*x2=32x1+24x2

     一天中因为检验员失误所造成的损失是：

     (8*25*0.02*x1+8*15*0.05*x2)*2=8x1+12x2

  约束条件：

    8*25*x1+8*15*x2>=1800

    8*25*x1<=1800

    8*15*x2<=1800

    X1>=0,x2>=0

即：

3 算法设计

li

粒子群初始化为二维数组，种群大小的行，2为列数，由所求的自变量为X1 和X2所得。

同样，粒子速度也为二维数组

粒子群算法执行过程中，加入约束条件检验函数入下：

代码：python语言

_date_ = '2019/5/27 20:16'
# -*- coding:utf-8 -*-
import numpy as np
import time
import matplotlib.pyplot as plt
import math
start = time.clock()
def getweight():
    # 惯性权重
    weight = 1
    return weight
def getlearningrate():
    # 分别是粒子的个体和社会的学习因子，也称为加速常数
    # lr = (0.1,1.2)
    lr = (0.49445, 1.49445)
    return lr
def getmaxgen():
    # 最大迭代次数
    maxgen = 20
    return maxgen
def getsizepop():
    # 种群规模
    sizepop = 20
    return sizepop
def getrangepop(i):
    # 粒子的位置的范围限制,x、y方向的限制相同
    # rangepop = (-2*math.pi , 2*math.pi)
    #rangepop = (-2,2)
    if (i ==0):
        rangepop = (0, 9)
    else:
        rangepop = (0, 15)
    return rangepop

def getrangespeed():
    # 粒子的速度范围限制
    # rangespeed = (-0.5,0.5)
    rangespeed = (-5, 5)
    return rangespeed

def func(x):
    # print("x的值和func")
    # print(x)
    if(issuitable(x)==0):
        y=10000
        return y
    # x输入粒子位置
    # y 粒子适应度值,求得最小函数值
    y =40*x[0]+36*x[1]
    # print(y)
    return y

#是否符合约束条件
def issuitable(x):
    # print("检验函数")
    if(5*x[0]+3*x[1]>=45):
        if((x[0]<=9.0)&(x[0]>=0.0)):
            if ((x[1] <=15.0)&( x[1] >= 0.0)):
                return 1
    return 0

def initpopvfit(sizepop):
    pop = np.zeros((sizepop,2))
    v = np.zeros((sizepop,2))
    fitness = np.zeros(sizepop)
    for i in range(sizepop):
        flag = 0
        while True:
            for j in range(2):
                pop[i][j] = (np.random.rand() - 0.5) * getrangepop(j)[1]*2
                pop[i][j] = round(pop[i][j])
                pop[pop < getrangepop(j)[0]] = getrangepop(j)[0]
                pop[pop > getrangepop(j)[1]] = getrangepop(j)[1]
                v[i][j] = np.random.rand()-0.5
                # 计算这个粒子的适应度
            flag = issuitable(pop[i])
            if flag == 1:
                fitness[i] = func(pop[i])
                break
    # print("初始fitness")
    # print(fitness)
    return pop,v,fitness

def getinitbest(fitness,pop):
    # 群体最优的粒子位置及其适应度值
    gbestpop,gbestfitness = pop[fitness.argmin()].copy(),fitness.min()
    #个体最优的粒子位置及其适应度值,使用copy()使得对pop的改变不影响pbestpop，pbestfitness类似
    pbestpop,pbestfitness = pop.copy(),fitness.copy()
    # print("全局最优")
    # print(gbestpop)
    # print(gbestfitness)
    return gbestpop,gbestfitness,pbestpop,pbestfitness

w = getweight()
# 获取惯性权重
lr = getlearningrate()
#分别是粒子的个体和社会的学习率
maxgen = getmaxgen()
# 获得最大迭代次数
sizepop = getsizepop()
# 获取种群规模
# rangepop = getrangepop()
# 粒子的位置的范围限制,x、y方向的限制相同
rangespeed = getrangespeed()
# 粒子的速度范围限制
pop,v,fitness = initpopvfit(sizepop)
#初始化种群
gbestpop,gbestfitness,pbestpop,pbestfitness = getinitbest(fitness,pop)
# 群体最优的粒子位置及其适应度值，个体最优的粒子位置及其适应度值
result = np.zeros(maxgen)
before = np.zeros(2)
beforev = np.zeros(2)
#记录粒子的自变量的值
#返回来一个maxgen大小的一维数组；
for u in range(maxgen):
    for i in range(sizepop):
        flag = 0
        for o in range(2):
            before[o] = pop[i][o]
            beforev=pop[i][o]
        while True:
            for j in range(2):
                v[i][j] = v[i][j]+lr[0] * np.random.rand() * (pbestpop[i][j] - pop[i][j]) + lr[1] * np.random.rand() * (gbestpop[j] - pop[i][j])
                v[v < rangespeed[0]] = rangespeed[0]
                v[v > rangespeed[1]] = rangespeed[1]
                # 粒子位置更新,需要检查粒子的约束条件是否符合
                pop[i][j] = v[i][j] + pop[i][j]
                pop[i][j]=round(pop[i][j])
                # int是向下取整
                pop[pop < getrangepop(j)[0]] = getrangepop(j)[0]
                pop[pop > getrangepop(j)[1]] = getrangepop(j)[1]
            flag = issuitable(pop[i])
            if flag == 1:
                break
            else:
                pop[i] = before
                v[i]=beforev
        #适应度更新
        for j in range(sizepop):
            fitness[j] = func(pop[j])
        for j in range(sizepop):
            if ((fitness[j] < pbestfitness[j])):
                pbestfitness[j] = fitness[j]
                pbestpop[j] = pop[j].copy()
            if(pbestfitness.min() < gbestfitness):
                gbestfitness = pbestfitness.min()
                gbestpop = pop[pbestfitness.argmin()].copy()
        result[u] = gbestfitness

print(result)
plt.plot(result)
plt.show()
end = time.clock()
print("time")
print((end-start)*1000.0)

由于自变量的取值空间较小，可用穷举法进行比较：

穷举法代码：

_date_ = '2019/5/28 7:43'
#利用穷举计算整数规划的值
import numpy as np
import time
start = time.clock()
num=0
max=0
m=0
n=0
b = np.zeros((160,3))
for x1 in range(10):
    for x2 in range(16):
        if (5 * x1 + 3 * x2>= 45):
            b[num][0]=x1
            b[num][1] = x2
            b[num][2] = 40*x1+36*x2
            num+=1
print("num")
print(num)
min=b[0][2]
m=b[0][0]
n=b[0][1]

for i in range(num):
    # print("b")
    # print(b)
    # print("b[num][2]")
    # print(b[0][2])
    if(b[i][2]<min):
        min=b[i][2]
        m=b[i][0]
        n=b[i][1]
        print(min)
print(min)
print(m)
print(n)
end = time.clock()
print("time")
print((end-start)*1000.0)
#ms为单位

对比结果：

用时竟然比穷举法还长好多，最后，尝试失败

可以找到使得粒子位置更符合约束条件的方法效果应该会好很多

出现的问题：

1 速度区间问题

   由用户设定用来限制粒子的速度

      区间太小加上约束条件容易陷入死循环

2 约束条件问题

约束条件太复杂的话，使用随机数更新的粒子位置很难满足约束条件