本文主要是介绍数学建模第二日,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!
今天换算法了!
`import geatpy as ga # 解遗传算法的库
from pymprog import * # 解运输问题的库
import numpy as np
import xlrd
-------------------解码函数---------------------------
def decode(Y, A, B, m):
# step1
List1 = np.multiply(np.array(Y), np.array(A)) # 计算 y*a
sum1 = 0 # 再把他们累加起来
for i in List1:sum1 += isum2 = 0
for j in B:sum2 += j# step2
while sum1 < sum2:k = np.random.randint(0, m)if Y[k] == 0:sum1 += A[k]Y[k] = 1else:k = np.random.randint(1, m + 1)continue
# else:# print("算法终止")return Y
-----------------------适值函数-------------------------------
固定成本函数
def Fixed_cost(List2):
fixed_cost = 0 # 再进行累加计算 计算建厂的固定成本
for i in List2:
fixed_cost += i
return fixed_cost
运输成本函数
def report(Y_decode, A, B): # 返回(产地,销地):运价;产地:产量; 销地:销量
# 把目标函数简化成产大于销的运输问题
a = (“A1”, “A2”, “A3”, “A4”, “A5”, “A6”, “A7”, “A8”, “A9”, “A10”) # 工厂(产地)
b = (“B1”, “B2”, “B3”, “B4”, “B5”, “B6”, “B7”, “B8”, “B9”, “B10”,
“B11”, “B12”, “B13”, “B14”, “B15”, “B16”, “B17”, “B18”, “B19”, “B20”, “B21”) # 需求地(销地)
datas = dict() # 定义一个空字典 (产地,销地):运价
datac = dict() # 定义一个空字典 产地:产量
datax = dict() # 定义一个空字典 销地:销量
List3 = np.multiply(np.array(Y_decode), np.array(A)) # 计算 y*a 产量
# (产地,销地):运价
Data = xlrd.open_workbook('单位运输费用.xlsx') # 读取Excel运价列表
table = Data.sheets()[0] # 读取Excel列表的sheet1
for i in range(0, 10): # 先从十行里面拿出一行来x = np.array(table.row_values(i))for j in range(0, 21): # 再从21列中取一列data = x[j]datas[a[i], b[j]] = data # (产地,销地):运价# 产地:产量
y = List3for i in range(0, 10):datac[a[i]] = y[i] # 产地:产量# 销地:销量
for j in range(0, 21):datax[b[j]] = B[j] # 销地:销量
return (datas, datac, datax)
def Transport_cost(datas, datac, datax):
a = (“A1”, “A2”, “A3”, “A4”, “A5”, “A6”, “A7”, “A8”, “A9”, “A10”) # 工厂(产地)
b = (“B1”, “B2”, “B3”, “B4”, “B5”, “B6”, “B7”, “B8”, “B9”, “B10”,
“B11”, “B12”, “B13”, “B14”, “B15”, “B16”, “B17”, “B18”, “B19”, “B20”, “B21”) # 需求地(销地)
# 开始模型计算
X = np.zeros((10, 21)) # 生成10行21列的0矩阵
cost = 0.0
begin(‘Transport’)
x = var(‘x’, datas.keys()) # 调运方案
minimize(sum(datas[i, j] * x[i, j] for (i, j) in datas.keys()), ‘Cost’) # 总运费最少
for i in datac.keys(): # 产地产量约束
sum(x[i, j] for j in datax.keys()) == datac[i]
for j in datax.keys(): # 销地销量约束
sum(x[i, j] for i in datac.keys()) == datax[j]
solve()
for (i, j) in datas.keys():# print('x[i, j].primal',x[i, j].primal)if x[i, j].primal > 0 and datas[i, j] != 0:# print("产地:%s -> 销地:%s 运输量:%-2d 运价:%2d" % (i, j, int(x[i, j].primal), int(datas[i, j])))X[a.index(i)][b.index(j)] = int(x[i, j].primal)
cost = int(vobj())
# print("X", X)
# print("总费用:%d" % int(vobj()))
end()
return (X, cost)
-----------------------主函数开始-------------------------------------
def main():
# 需求点的需求量
B = [6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 0]
# 工厂容量(生产能力)
A = [30, 40, 50, 55, 45, 48, 36, 50, 60, 35]# 各待选工厂的固定费用
D = [65, 50, 40, 30, 20, 50, 45, 30, 35, 25]# 从待选工厂到需求点的单位运输费用
Data = xlrd.open_workbook('单位运输费用.xlsx')
table = Data.sheets()[0]Y = [] # 定义变量码
Y_decode = [] # 定义解码生成的变量码
m = 10 # 工厂个数
T = 100 # 迭代次数
decode_after_chrom2 = np.zeros((60, 10)) # 存储60条Y_decode
datas = dict() # 定义一个空字典 (产地,销地):运价
datac = dict() # 定义一个空字典 产地:产量
datax = dict() # 定义一个空字典 销地:销量
fixed_cost = 0.0 # 初始固定成本
X = np.zeros((10, 21)) # 运量矩阵
transport_cost = 0.0 # 初始运输费用
cost = 0.0 # 初始总费用best_y = np.zeros((1, 10)) # 所有代中的最优选址方案
min_t = 0 # 最优的方案出现的代数
all_min_cost = 1000 # 所有代中最小总费用
all_min_X = np.zeros((10, 21)) # 所有代中最小总费用的运量矩阵
FitnV = np.zeros((1, 60)) # 适应度列表
NewChrlx = [] # 轮盘赌存储位置
New_chrom1 = np.zeros((60, 10))
New_chrom2 = np.zeros((60, 10))# ---------------------遗传编码---------------------------# 生成附加码
chrom1 = ga.crtpp(60, 10, 10) # 种群大小,染色体长度,种群最大元素
print(chrom1)# 生成变量码
FieldDR = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0],[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]])
chrom2 = ga.crtip(60, FieldDR) # 种群大小,上界下界区域描述器
print(chrom2)
add_code = np.zeros((1, 10)) # 定义父代单行附加码
new_add_code = np.zeros((1, 10)) # 定义子代单行附加码
variable_code = np.zeros((1, 10)) # 定义父代单行变量码
new_variable_code = np.zeros((1, 10)) # 定义子代单行变量码
index_v = 0 # 定义初始索引值
begin_New_chrom2 = np.zeros((60, 10)) # 定义经过交叉后的新的变量码种群
again_New_chrom1 = np.zeros((60, 10)) # 经过变异后的附加码种群
# tem_list = []
a = 0
for t in np.arange(0, T):min_cost = 1000 # 初始最小总费用min_X = np.zeros((10, 21)) # 初始最小总费用的运量矩阵min_i = 400 # 每一代中最小的运价在种群对应的位置better_y = np.zeros((1, 10)) # 每一代中的最有选址方案print("==========================第%d代种群======================" % t)for i in np.arange(0, 60):# print("第%d个染色体" % i)Y = chrom2[i]# print("-------------%d-----------------" % i)# print("Y",Y)Y_decode = np.array(decode(Y, A, B, m))decode_after_chrom2[i] = Y_decodeList2 = np.multiply(np.array(D), np.array(Y_decode)) # 计算 d*yfixed_cost = Fixed_cost(List2)print("fixed_cost", fixed_cost)(datas, datac, datax) = report(Y_decode, A, B)(X, transport_cost) = Transport_cost(datas, datac, datax)cost = fixed_cost + transport_costprint('第%d个染色体的cost' % i, cost)FitnV[0][i] = 1 / costif cost < min_cost:min_cost = costmin_X = Xmin_i = ibetter_y = Y_decode# print("第%d代种群中最小的cost"%t, min_cost)# print("最优运输方案", min_X)# print('对应的位置', min_i)NewChrlx = ga.rws(np.array(FitnV).T, 59) # 轮盘赌# print(NewChrlx)for i in NewChrlx:New_chrom1[i] = chrom1[i]New_chrom2[i] = decode_after_chrom2[i]for j in range(0, 60): # 把最大的适值得到染色体拿到新的种群中if j not in NewChrlx:New_chrom1[j] = chrom1[j]New_chrom2[j] = decode_after_chrom2[j]# 遗传操作Pc = 0.6 # 交叉概率Pm = 0.1 # 变异概率New_chrom1 = ga.xovpm(New_chrom1, 0.6)# print("交叉后的附加码New_chrom1",New_chrom1)# print("chrom1",chrom1)for change_1 in np.arange(0, 60):# print('------------修改chrom2第%d次-----------'%change_1)new_add_code = New_chrom1[change_1]add_code = chrom1[change_1]# print(type(add_code))variable_code = New_chrom2[change_1]for change_2 in np.arange(0, 10):index_v = add_code.tolist().index(new_add_code[change_2])# print("index_v",index_v)new_variable_code[0, change_2] = variable_code[index_v]begin_New_chrom2[change_1] = new_variable_code# print("begin_New_chrom2[change_1]",begin_New_chrom2[change_1])new_variable_code = np.zeros((1, 10))# print("交叉后的变量码begin_New_chrom2",begin_New_chrom2)# 变异操作ret = np.random.random()if ret > Pm:for mutation_v in np.arange(0, 60):new_add_code = New_chrom1[mutation_v]# print('未翻转的附加码new_add_code',new_add_code)random_1 = np.random.randint(0, 10)random_2 = np.random.randint(0, 10)a = np.abs(random_1 - random_2)tem_list = np.zeros(a + 1)if a != 0:if random_2 > random_1:for i, k in zip(np.arange(0, np.abs(random_2 - random_1) + 1), np.arange(random_1, random_2 + 1)):tem_list[i] = new_add_code[k]# print("未翻转的单行tem_list",tem_list)# print("未翻转的tem_list类型",type(tem_list))tem_list = tem_list.tolist()# print("转换为列表的tem_list", tem_list)tem_list.reverse()# print("翻转后的tem_list",tem_list)for i, k in zip(np.arange(0, np.abs(random_2 - random_1) + 1).tolist(),np.arange(random_1, random_2 + 1).tolist()):new_add_code[k] = tem_list[i]# print("翻转后的附加码",new_add_code)again_New_chrom1[mutation_v] = new_add_codeelse:for i, k in zip(np.arange(0, np.abs(random_2 - random_1) + 1), np.arange(random_2, random_1 + 1)):tem_list[i] = new_add_code[k]# print("未翻转的单行tem_list",tem_list)# print("未翻转的tem_list类型",type(tem_list))tem_list = tem_list.tolist()# print("转换为列表的tem_list",tem_list)tem_list.reverse()# print("tem_list",tem_list)for i, k in zip(np.arange(0, np.abs(random_2 - random_1) + 1).tolist(),np.arange(random_2, random_1 + 1).tolist()):new_add_code[k] = tem_list[i]again_New_chrom1[mutation_v] = new_add_code# print("循环中的again_New_chrom1", again_New_chrom1)else:again_New_chrom1[mutation_v] = new_add_code# print("else循环中的again_New_chrom1", again_New_chrom1)# print("New_chrom1",New_chrom1)# print('变异后的附加码again_New_chrom1',again_New_chrom1)chrom1 = again_New_chrom1chrom2 = begin_New_chrom2print('遗传操作后的chrom2', chrom2)print("第%d代种群中最小的cost" % t, min_cost)print("最优运输方案", min_X)print('对应的位置', min_i)if all_min_cost > min_cost:all_min_cost = min_costall_min_X = min_Xbest_y = better_ymin_t = t
# New_chrom1[60]= chrom1[min_i]
# New_chrom2[60] = decode_after_chrom2[min_i]
# print("FitnV",np.array(FitnV).T.shape)
print("所有代中种群中最小的cost", all_min_cost)
print("所有代中的最优运输方案", all_min_X)
print("最优的选址方案", best_y)
print("最优的选址方案出现在第%d代" % min_t)
if __name__ == "__main__":
main()
`
import numpy as np
import matplotlib.pyplot as plt
import pdb"旅行商问题 ( TSP , Traveling Salesman Problem )"
coordinates = np.array([[565.0,575.0],[25.0,185.0],[345.0,750.0],[945.0,685.0],[845.0,655.0],[880.0,660.0],[25.0,230.0],[525.0,1000.0],[580.0,1175.0],[650.0,1130.0],[1605.0,620.0],[1220.0,580.0],[1465.0,200.0],[1530.0, 5.0],[845.0,680.0],[725.0,370.0],[145.0,665.0],[415.0,635.0],[510.0,875.0],[560.0,365.0],[300.0,465.0],[520.0,585.0],[480.0,415.0],[835.0,625.0],[975.0,580.0],[1215.0,245.0],[1320.0,315.0],[1250.0,400.0],[660.0,180.0],[410.0,250.0],[420.0,555.0],[575.0,665.0],[1150.0,1160.0],[700.0,580.0],[685.0,595.0],[685.0,610.0],[770.0,610.0],[795.0,645.0],[720.0,635.0],[760.0,650.0],[475.0,960.0],[95.0,260.0],[875.0,920.0],[700.0,500.0],[555.0,815.0],[830.0,485.0],[1170.0, 65.0],[830.0,610.0],[605.0,625.0],[595.0,360.0],[1340.0,725.0],[1740.0,245.0]])#得到距离矩阵的函数
def getdistmat(coordinates):num = coordinates.shape[0] #52个坐标点distmat = np.zeros((52,52)) #52X52距离矩阵for i in range(num):for j in range(i,num):distmat[i][j] = distmat[j][i]=np.linalg.norm(coordinates[i]-coordinates[j])return distmatdef initpara():alpha = 0.99t = (1,100)markovlen = 1000return alpha,t,markovlen
num = coordinates.shape[0]
distmat = getdistmat(coordinates) #得到距离矩阵solutionnew = np.arange(num)
#valuenew = np.max(num)solutioncurrent = solutionnew.copy()
valuecurrent =99000 #np.max这样的源代码可能同样是因为版本问题被当做函数不能正确使用,应取一个较大值作为初始值
#print(valuecurrent)solutionbest = solutionnew.copy()
valuebest = 99000 #np.maxalpha,t2,markovlen = initpara()
t = t2[1]result = [] #记录迭代过程中的最优解
while t > t2[0]:for i in np.arange(markovlen):#下面的两交换和三角换是两种扰动方式,用于产生新解if np.random.rand() > 0.5:# 交换路径中的这2个节点的顺序# np.random.rand()产生[0, 1)区间的均匀随机数while True:#产生两个不同的随机数loc1 = np.int(np.ceil(np.random.rand()*(num-1)))loc2 = np.int(np.ceil(np.random.rand()*(num-1)))## print(loc1,loc2)if loc1 != loc2:breaksolutionnew[loc1],solutionnew[loc2] = solutionnew[loc2],solutionnew[loc1]else: #三交换while True:loc1 = np.int(np.ceil(np.random.rand()*(num-1)))loc2 = np.int(np.ceil(np.random.rand()*(num-1))) loc3 = np.int(np.ceil(np.random.rand()*(num-1)))if((loc1 != loc2)&(loc2 != loc3)&(loc1 != loc3)):break# 下面的三个判断语句使得loc1<loc2<loc3if loc1 > loc2:loc1,loc2 = loc2,loc1if loc2 > loc3:loc2,loc3 = loc3,loc2if loc1 > loc2:loc1,loc2 = loc2,loc1#下面的三行代码将[loc1,loc2)区间的数据插入到loc3之后tmplist = solutionnew[loc1:loc2].copy()solutionnew[loc1:loc3-loc2+1+loc1] = solutionnew[loc2:loc3+1].copy()solutionnew[loc3-loc2+1+loc1:loc3+1] = tmplist.copy() valuenew = 0for i in range(num-1):valuenew += distmat[solutionnew[i]][solutionnew[i+1]]valuenew += distmat[solutionnew[0]][solutionnew[51]]# print (valuenew)if valuenew<valuecurrent: #接受该解#更新solutioncurrent 和solutionbestvaluecurrent = valuenewsolutioncurrent = solutionnew.copy()if valuenew < valuebest:valuebest = valuenewsolutionbest = solutionnew.copy()else:#按一定的概率接受该解if np.random.rand() < np.exp(-(valuenew-valuecurrent)/t):valuecurrent = valuenewsolutioncurrent = solutionnew.copy()else:solutionnew = solutioncurrent.copy()t = alpha*tresult.append(valuebest)#print (t) #程序运行时间较长,打印t来监视程序进展速度print(valuebest)
#用来显示结果
plt.plot(np.array(result))
plt.ylabel("bestvalue")
plt.xlabel("t")
plt.show()
https://blog.csdn.net/weixin_44086712/article/details/107413339
参考以上三点!
这篇关于数学建模第二日的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!