diff --git a/.gitignore b/.gitignore index a5a8fb6..02618ef 100644 --- a/.gitignore +++ b/.gitignore @@ -8,4 +8,5 @@ Scripts/ *.cfg *.pkl *.pth -*.m \ No newline at end of file +*.m +opt/ \ No newline at end of file diff --git a/base_optimizer/optimizer_common.py b/base_optimizer/optimizer_common.py index a6c7232..365d67c 100644 --- a/base_optimizer/optimizer_common.py +++ b/base_optimizer/optimizer_common.py @@ -22,6 +22,7 @@ import pandas as pd import matplotlib.pyplot as plt import matplotlib import traceback +import openpyxl matplotlib.use('TkAgg') @@ -67,7 +68,16 @@ t_fix_camera_check = 0.12 # 固定相机检测时间 T_pp, T_tr, T_nc, T_pl = 2, 5, 25, 0 # 时间参数 （数据拟合获得） -Fit_cy, Fit_nz, Fit_pu, Fit_pl, Fit_mv = 0.326, 0.870, 0.159, 0.041, 0.001 +# Fit_cy, Fit_nz, Fit_pu, Fit_pl, Fit_mv = 0.326, 0.870, 0.159, 0.041, 0.030 +Fit_cy, Fit_nz, Fit_pu, Fit_pl, Fit_mv = 0.326, 0.870, 0.159, 0.041, 0.035 + + +class Point: + def __init__(self, _x, _y, _r=0, _h=None): + self.x = _x + self.y = _y + self.r = _r + self.h = _h class OptResult: @@ -108,9 +118,12 @@ class OptInfo: print(f'-Pick time: {self.pickup_time: .3f}, Pick distance: {self.pickup_distance: .3f}') print(f'-Place time: {self.place_time: .3f}, Place distance: {self.place_distance: .3f}') print( - f'-Round time: {self.total_time - self.place_time - self.place_time: .3f}, Round distance: ' + f'-Round time: {self.total_time - self.operation_time - self.pickup_time - self.place_time: .3f}, Round distance: ' f'{self.total_distance - self.pickup_distance - self.place_distance: .3f}') + print(f'-Round & place time per cycle: {(self.total_time - self.pickup_time - self.operation_time) * 1000.0 / (self.cycle_counter + 1e-10): .3f}, ', end='') + print(f'-Round & place distance per cycle: {(self.total_distance - self.pickup_distance) / (self.cycle_counter + 1e-10): .3f}') + minutes, seconds = int(self.total_time // 60), int(self.total_time) % 60 millisecond = int((self.total_time - minutes * 60 - seconds) * 60) @@ -856,3 +869,10 @@ def random_division(num, div): def list_range(start, end=None): return list(range(start)) if end is None else list(range(start, end)) + + +def kth_indices_partition(num, kth): + if len(num) > kth: + return np.argpartition(num, kth) + else: + return np.array(range(len(num))) \ No newline at end of file diff --git a/base_optimizer/optimizer_interface.py b/base_optimizer/optimizer_interface.py index a524dec..a198801 100644 --- a/base_optimizer/optimizer_interface.py +++ b/base_optimizer/optimizer_interface.py @@ -15,26 +15,22 @@ def base_optimizer(machine_index, pcb_data, component_data, feeder_data, params, if params.machine_optimizer == 'cell-division': # 基于元胞分裂的遗传算法 component_result, cycle_result, feeder_slot_result = optimizer_celldivision(pcb_data, component_data) - placement_result, head_sequence = greedy_placement_route_generation(component_data, pcb_data, component_result, - cycle_result, feeder_slot_result) + placement_result, head_sequence = place_allocate_sequence_route_generation(component_data, pcb_data, + component_result, cycle_result, + feeder_slot_result) elif params.machine_optimizer == 'feeder-priority': # 基于基座扫描的供料器优先算法 component_result, cycle_result, feeder_slot_result = feeder_priority_assignment(component_data, pcb_data, feeder_data) - placement_result, head_sequence = greedy_placement_route_generation(component_data, pcb_data, component_result, - cycle_result, feeder_slot_result) - # placement_result, head_sequence = beam_search_route_generation(component_data, pcb_data, component_result, - # cycle_result, feeder_slot_result) - # placement_result, head_sequence = scan_based_placement_route_generation(component_data, pcb_data, - # component_result, cycle_result, - # feeder_slot_result) - + placement_result, head_sequence = scan_based_placement_route_generation(component_data, pcb_data, + component_result, cycle_result, + feeder_slot_result) elif params.machine_optimizer == 'hybrid-genetic': # 基于拾取组的混合遗传算法 component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_hybrid_genetic( pcb_data, component_data, hinter=hinter) elif params.machine_optimizer == 'aggregation': # 基于batch-level的整数规划 + 启发式算法 component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_aggregation( - component_data, pcb_data) + component_data, pcb_data, hinter=hinter) elif params.machine_optimizer == 'genetic-scanning': component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_genetic_scanning( component_data, pcb_data, hinter=hinter) @@ -43,15 +39,16 @@ def base_optimizer(machine_index, pcb_data, component_data, feeder_data, params, component_data, pcb_data, hinter=hinter) elif params.machine_optimizer == "two-phase": component_result, feeder_slot_result, cycle_result = gurobi_optimizer(pcb_data, component_data, feeder_data, - initial=True, partition=True, + initial=True, partition=False, reduction=True, hinter=hinter) - placement_result, head_sequence = greedy_placement_route_generation(component_data, pcb_data, component_result, - cycle_result, feeder_slot_result) + placement_result, head_sequence = place_allocate_sequence_route_generation(component_data, pcb_data, + component_result, cycle_result, + feeder_slot_result) else: raise 'machine optimizer method ' + params.method + ' is not existed' - print('----- Placement machine ' + str(machine_index) + ' ----- ') + # print('----- Placement machine ' + str(machine_index) + ' ----- ') opt_res = OptResult(component_result, cycle_result, feeder_slot_result, placement_result, head_sequence) # 估算贴装用时 info = placement_info_evaluation(component_data, pcb_data, opt_res, hinter=False) @@ -67,4 +64,6 @@ def base_optimizer(machine_index, pcb_data, component_data, feeder_data, params, f'result/{params.filename[:-4]}-{params.line_optimizer}-M0{machine_index} {params.save_suffix}', component_data, pcb_data, opt_res) + # output_optimize_result(f'{params.filename[:-4]}', component_data, pcb_data, opt_res) + return info diff --git a/base_optimizer/result_analysis.py b/base_optimizer/result_analysis.py index 22431b5..9187a04 100644 --- a/base_optimizer/result_analysis.py +++ b/base_optimizer/result_analysis.py @@ -355,15 +355,12 @@ def output_optimize_result(file_path, component_data, pcb_data, optimizer_result if 'desc' not in output_data.columns: column_index = int(np.where(output_data.columns.values.reshape(-1) == 'part')[0][0]) output_data.insert(loc=column_index + 1, column='desc', value='') - file_dir = file_path[:file_path.rfind('/') + 1] - if not os.path.exists(file_dir): - os.makedirs(file_dir) - output_data.to_excel(file_path + '.xlsx', sheet_name='tb1', float_format='%.3f', na_rep='') + output_data.to_csv('result/' + file_path + '.txt', sep='\t', float_format='%.3f', header=False, index=False) def optimization_assign_result(component_data, pcb_data, optimizer_result, nozzle_hinter=False, component_hinter=False, - feeder_hinter=False): + feeder_hinter=False, placement_hinter=False): if nozzle_hinter: columns = ['H{}'.format(i + 1) for i in range(max_head_index)] + ['cycle'] @@ -428,6 +425,29 @@ def optimization_assign_result(component_data, pcb_data, optimizer_result, nozzl print(feedr_assign) print('') + if placement_hinter: + columns = ['H{}'.format(i + 1) for i in range(max_head_index)] + ['cycle'] + + placement_assign = pd.DataFrame(columns=columns) + for cycle, _ in enumerate(optimizer_result.placement_assign): + placement_assign.loc[cycle, 'cycle'] = 1 + for head in range(max_head_index): + point = optimizer_result.placement_assign[cycle][head] + if point != -1: + placement_assign.loc[cycle, 'H{}'.format(head + 1)] = 'P{}'.format(point) + else: + placement_assign.loc[cycle, 'H{}'.format(head + 1)] = '' + + headseq_assign = pd.DataFrame(columns=columns) + for cycle, headseq in enumerate(optimizer_result.head_sequence): + headseq_assign.loc[cycle, 'cycle'] = 1 + for head in range(len(headseq)): + headseq_assign.loc[cycle, 'H{}'.format(head + 1)] = 'H{}'.format(headseq[head]) + + print(placement_assign) + print(headseq_assign) + print('') + def placement_info_evaluation(component_data, pcb_data, optimizer_result, hinter=False): # === 优化结果参数 === @@ -555,28 +575,43 @@ def placement_info_evaluation(component_data, pcb_data, optimizer_result, hinter # 贴装路径 if optimizer_result.placement_assign and optimizer_result.head_sequence: + head_angle = [0 for _ in range(max_head_index)] for head in optimizer_result.head_sequence[cycle]: index = optimizer_result.placement_assign[cycle][head] if index == -1: continue mount_pos.append([pcb_data.iloc[index]['x'] - head * head_interval + stopper_pos[0], pcb_data.iloc[index]['y'] + stopper_pos[1]]) - mount_angle.append(pcb_data.iloc[index]['r']) + head_angle[head] = pcb_data.iloc[index]['r'] # 单独计算贴装路径 for cntPoints in range(len(mount_pos) - 1): info.place_distance += max(abs(mount_pos[cntPoints][0] - mount_pos[cntPoints + 1][0]), abs(mount_pos[cntPoints][1] - mount_pos[cntPoints + 1][1])) + if mount_pos[0][0] < mount_pos[-1][0]: + mount_pos = reversed(mount_pos) + # 考虑R轴预旋转，补偿同轴角度转动带来的额外贴装用时 - info.operation_time += head_rotary_time(mount_angle[0]) # 补偿角度转动带来的额外贴装用时 info.operation_time += t_nozzle_put * nozzle_put_counter + t_nozzle_pick * nozzle_pick_counter for idx, pos in enumerate(mount_pos): info.operation_time += t_place - move_time = max(axis_moving_time(cur_pos[0] - pos[0], 0), axis_moving_time(cur_pos[1] - pos[1], 1)) + if idx == 0: + move_time = max(axis_moving_time(cur_pos[0] - pos[0], 0), + axis_moving_time(cur_pos[1] - pos[1], 1)) info.round_time += move_time else: + + cur_head = optimizer_result.head_sequence[cycle][idx] + side_head = cur_head - 1 if cur_head % 2 else cur_head + 1 + if optimizer_result.head_sequence[cycle][idx - 1] != side_head: + move_time = max(axis_moving_time(cur_pos[0] - pos[0], 0), + axis_moving_time(cur_pos[1] - pos[1], 1)) + else: + move_time = max(axis_moving_time(cur_pos[0] - pos[0], 0), + axis_moving_time(cur_pos[1] - pos[1], 1), + head_rotary_time(head_angle[cur_head] - head_angle[side_head])) info.place_time += move_time info.total_distance += max(abs(cur_pos[0] - pos[0]), abs(cur_pos[1] - pos[1])) @@ -594,5 +629,3 @@ def placement_info_evaluation(component_data, pcb_data, optimizer_result, hinter return info - - diff --git a/base_optimizer/smopt_aggregation.py b/base_optimizer/smopt_aggregation.py index 58385d4..8a687e9 100644 --- a/base_optimizer/smopt_aggregation.py +++ b/base_optimizer/smopt_aggregation.py @@ -1,12 +1,12 @@ from base_optimizer.optimizer_common import * - +from base_optimizer.smtopt_route import * def list_range(start, end=None): return list(range(start)) if end is None else list(range(start, end)) @timer_wrapper -def optimizer_aggregation(component_data, pcb_data): +def optimizer_aggregation(component_data, pcb_data, hinter=True): # === phase 0: data preparation === M = 1000 # a sufficient large number a, b = 1, 6 # coefficient @@ -43,7 +43,7 @@ def optimizer_aggregation(component_data, pcb_data): # === phase 1: mathematical model solver === mdl = Model('SMT') - mdl.setParam('OutputFlag', 0) + mdl.setParam('OutputFlag', hinter) # === Decision Variables === # the number of components of type i that are placed by nozzle type j on placement head k @@ -104,8 +104,7 @@ def optimizer_aggregation(component_data, pcb_data): mdl.setParam("TimeLimit", 100) mdl.optimize() - - if mdl.Status == GRB.OPTIMAL: + if mdl.Status == GRB.OPTIMAL or mdl.Status == GRB.TIME_LIMIT: print('total cost = {}'.format(mdl.objval)) # convert cp model solution to standard output @@ -160,25 +159,9 @@ def optimizer_aggregation(component_data, pcb_data): feeder_slot_result = feeder_assignment(component_data, pcb_data, component_result, cycle_result) # === phase 2: heuristic method === - mount_point_pos = defaultdict(list) - for pcb_idx, data in pcb_data.iterrows(): - part = data['part'] - part_index = component_data[component_data['part'] == part].index.tolist()[0] - mount_point_pos[part_index].append([data['x'], data['y'], pcb_idx]) - - for index_ in mount_point_pos.keys(): - mount_point_pos[index_].sort(key=lambda x: (x[1], x[0])) - - for cycle_idx, _ in enumerate(cycle_result): - for _ in range(cycle_result[cycle_idx]): - placement_result.append([-1 for _ in range(max_head_index)]) - for head in range(max_head_index): - if component_result[cycle_idx][head] == -1: - continue - index_ = component_result[cycle_idx][head] - placement_result[-1][head] = mount_point_pos[index_][-1][2] - mount_point_pos[index_].pop() - head_sequence.append(dynamic_programming_cycle_path(pcb_data, placement_result[-1], feeder_slot_result[cycle_idx])) + placement_result, head_sequence = greedy_level_placing_route_generation(component_data, pcb_data, + component_result, cycle_result, + feeder_slot_result) else: warnings.warn('No solution found!', UserWarning) diff --git a/base_optimizer/smopt_celldivision.py b/base_optimizer/smopt_celldivision.py index 789f2fc..c6851cc 100644 --- a/base_optimizer/smopt_celldivision.py +++ b/base_optimizer/smopt_celldivision.py @@ -33,7 +33,7 @@ def convert_cell_2_result(pcb_data, component_data, component_cell, population): wl = [0 for _ in range(max_head_index)] # workload - e1, e2, e3 = 1, 0.5, 1. / 6 + e1, e2, e3 = 1, 2, 1. / 6 component_result, cycle_result, feeder_slot_result = [], [], [] for index in population: @@ -101,12 +101,16 @@ def optimizer_celldivision(pcb_data, component_data, hinter=True): golden_section = 0.618 # 获取元件元胞 - point_num = len(pcb_data) - component_cell = pd.DataFrame({'index': np.arange(len(component_data)), 'points': np.zeros(len(component_data), dtype=int)}) - for point_cnt in range(point_num): - part = pcb_data.loc[point_cnt, 'part'] - index = np.where(component_data['part'].values == part) - component_cell.loc[index[0], 'points'] += 1 + feeder_num = sum(component_data['fdn']) + component_cell = pd.DataFrame({'index': np.arange(feeder_num), 'points': np.zeros(feeder_num, dtype=int)}) + cell_index = 0 + for part_index, data in component_data.iterrows(): + total_points, division_points = data.points, math.ceil(data.points / data.fdn) + for _ in range(data.fdn): + component_cell.loc[cell_index, 'index'] = part_index + component_cell.loc[cell_index, 'points'] = min(division_points, total_points) + total_points -= division_points + cell_index += 1 component_cell = component_cell[~component_cell['points'].isin([0])] # component_cell.sort_values(by = "points" , inplace = True, ascending = False) diff --git a/base_optimizer/smopt_hybridgenetic.py b/base_optimizer/smopt_hybridgenetic.py index 04281ff..9e48157 100644 --- a/base_optimizer/smopt_hybridgenetic.py +++ b/base_optimizer/smopt_hybridgenetic.py @@ -1,5 +1,5 @@ from base_optimizer.optimizer_common import * - +from base_optimizer.smtopt_route import * def dynamic_programming_cycle_path(cycle_placement, cycle_points): head_sequence = [] @@ -231,8 +231,8 @@ def cal_individual_val(component_nozzle, component_point_pos, designated_nozzle, def convert_individual_2_result(component_data, component_point_pos, designated_nozzle, pickup_group, pickup_group_cycle, pair_group, feeder_lane, individual): + component_result, cycle_result, feeder_slot_result = [], [], [] - placement_result, head_sequence_result = [], [] # === 记录不同元件对应的槽位 === feeder_part_arrange = defaultdict(list) @@ -273,26 +273,7 @@ def convert_individual_2_result(component_data, component_point_pos, designated_ pickup_cycle_result[idx][head] -= cycle - component_point_index = defaultdict(int) - for cycle_set in range(len(cycle_result)): - for cycle in range(cycle_result[cycle_set]): - placement_result.append([-1 for _ in range(max_head_index)]) - mount_point = [[0, 0] for _ in range(max_head_index)] - for head in range(max_head_index): - part_index = component_result[cycle_set][head] - if part_index == -1: - continue - - part = component_data.iloc[part_index]['part'] - point_info = component_point_pos[part][component_point_index[part]] - - placement_result[-1][head] = point_info[2] - mount_point[head] = point_info[0:2] - - component_point_index[part] += 1 - head_sequence_result.append(dynamic_programming_cycle_path(placement_result[-1], mount_point)) - - return component_result, cycle_result, feeder_slot_result, placement_result, head_sequence_result + return component_result, cycle_result, feeder_slot_result @timer_wrapper @@ -510,11 +491,6 @@ def optimizer_hybrid_genetic(pcb_data, component_data, hinter=True): pop_val.append(val) # val is related to assembly time for _ in range(n_generations): - # idx = np.argmin(pop_val) - # if len(best_pop_val) == 0 or pop_val[idx] < best_pop_val[-1]: - # best_individual = copy.deepcopy(population[idx]) - # best_pop_val.append(pop_val[idx]) - # min-max convert max_val = 1.5 * max(pop_val) convert_pop_val = list(map(lambda v: max_val - v, pop_val)) @@ -565,6 +541,15 @@ def optimizer_hybrid_genetic(pcb_data, component_data, hinter=True): pbar.update(1) best_individual = population[np.argmin(pop_val)] + component_result, cycle_result, feeder_slot_result = convert_individual_2_result(component_data, + component_point_pos, + designated_nozzle, pickup_group, + pickup_group_cycle, pair_group, + feeder_lane, best_individual) + placement_result, head_sequence_result = place_cluster_greedy_route_generation(component_data, pcb_data, + component_result, cycle_result, + feeder_slot_result) + + + return component_result, cycle_result, feeder_slot_result, placement_result, head_sequence_result - return convert_individual_2_result(component_data, component_point_pos, designated_nozzle, pickup_group, - pickup_group_cycle, pair_group, feeder_lane, best_individual) diff --git a/base_optimizer/smopt_mathmodel.py b/base_optimizer/smopt_mathmodel.py index fbac9e3..f26663f 100644 --- a/base_optimizer/smopt_mathmodel.py +++ b/base_optimizer/smopt_mathmodel.py @@ -1,4 +1,5 @@ from base_optimizer.optimizer_common import * +from base_optimizer.smtopt_route import * def head_task_model(component_data, pcb_data, hinter=True): @@ -6,11 +7,10 @@ def head_task_model(component_data, pcb_data, hinter=True): mdl = Model('pick_route') mdl.setParam('Seed', 0) mdl.setParam('OutputFlag', hinter) # set whether output the debug information - mdl.setParam('TimeLimit', 600) + mdl.setParam('TimeLimit', 1000) H = max_head_index I = len(component_data) - S = len(component_data) K = len(pcb_data) nozzle_type, component_type = [], [] @@ -29,12 +29,16 @@ def head_task_model(component_data, pcb_data, hinter=True): M = 10000 CompOfNozzle = [[0 for _ in range(J)] for _ in range(I)] # Compatibility - component_point = [0 for _ in range(I)] + component_point, component_fdn = [0 for _ in range(I)], [0 for _ in range(I)] + for _, data in pcb_data.iterrows(): idx = component_data[component_data.part == data.part].index.tolist()[0] nozzle = component_data.iloc[idx].nz CompOfNozzle[idx][nozzle_type.index(nozzle)] = 1 component_point[idx] += 1 + component_fdn[idx] = component_data.iloc[idx].fdn + + S = sum(component_fdn) # objective related g = mdl.addVars(list_range(K), vtype=GRB.BINARY) @@ -94,7 +98,7 @@ def head_task_model(component_data, pcb_data, hinter=True): range(-(H - 1) * r, S) for k in range(K)) # feeder related - mdl.addConstrs(quicksum(f[s, i] for s in range(S)) <= 1 for i in range(I)) + mdl.addConstrs(quicksum(f[s, i] for s in range(S)) <= component_fdn[i] for i in range(I)) mdl.addConstrs(quicksum(f[s, i] for i in range(I)) <= 1 for s in range(S)) mdl.addConstrs( quicksum(x[i, k, h] * y[s, k, h] for h in range(H) for k in range(K)) >= f[s, i] for i in range(I) for s in range(S)) @@ -110,7 +114,7 @@ def head_task_model(component_data, pcb_data, hinter=True): # objective mdl.setObjective(Fit_cy * quicksum(g[k] for k in range(K)) + Fit_nz * quicksum( d[k, h] for h in range(H) for k in range(K)) + Fit_pu * quicksum( - e[s, k] for s in range(-(H - 1) * r, S) for k in range(K)) + Fit_mv * head_interval * quicksum(u[k] for k in range(K)), + e[s, k] for s in range(-(H - 1) * r, S) for k in range(K)) + Fit_mv * quicksum(u[k] for k in range(K)), GRB.MINIMIZE) mdl.optimize() @@ -119,19 +123,21 @@ def head_task_model(component_data, pcb_data, hinter=True): for k in range(K): if abs(g[k].x) < 1e-6: continue + component_assign, feeder_slot_assign = [-1 for _ in range(H)], [-1 for _ in range(H)] - component_result.append([-1 for _ in range(H)]) - feeder_slot_result.append([-1 for _ in range(H)]) - cycle_result.append(1) for h in range(H): for i in range(I): if abs(x[i, k, h].x) > 1e-6: - component_result[-1][h] = i + component_assign[h] = i for s in range(S): if abs(y[s, k, h].x) > 1e-6: - feeder_slot_result[-1][h] = slot_start + s * interval_ratio - 1 + feeder_slot_assign[h] = slot_start + s * interval_ratio - 1 + if sum(component_assign) != -H: + component_result.append(component_assign) + feeder_slot_result.append(feeder_slot_assign) + cycle_result.append(1) if hinter: print(component_result) @@ -143,7 +149,7 @@ def place_route_model(component_data, pcb_data, component_result, feeder_slot_re mdl = Model('place_route') mdl.setParam('Seed', 0) mdl.setParam('OutputFlag', hinter) # set whether output the debug information - mdl.setParam('TimeLimit', 10) + mdl.setParam('TimeLimit', 1000) component_type = [] for _, data in component_data.iterrows(): @@ -227,7 +233,7 @@ def place_route_model(component_data, pcb_data, component_result, feeder_slot_re else: # there are components on the head mdl.addConstrs(quicksum(w[p, q, k, a] for a in A_from(h) for q in range(P)) + quicksum( - w[q, p, k, a] for a in A_to(h) for q in range(P)) + y[p, k, h] + z[p, k, h] <= 4 * + w[q, p, k, a] for a in A_to(h) for q in range(P)) + y[p, k, h] + z[p, k, h] <= 2 * CompOfPoint[component_result[k][h]][p] for p in range(P)) # each head corresponds to a maximum of one point in each cycle @@ -362,7 +368,9 @@ def place_route_model(component_data, pcb_data, component_result, feeder_slot_re def optimizer_mathmodel(component_data, pcb_data, hinter=True): component_result, cycle_result, feeder_slot_result = head_task_model(component_data, pcb_data, hinter) - placement_result, head_sequence = place_route_model(component_data, pcb_data, component_result, feeder_slot_result) - # placement_result, head_sequence = greedy_placement_route_generation(component_data, pcb_data, component_result, - # cycle_result) + # placement_result, head_sequence = place_route_model(component_data, pcb_data, component_result, feeder_slot_result) + placement_result, head_sequence = place_allocate_sequence_route_generation(component_data, pcb_data, + component_result, cycle_result, + feeder_slot_result) + return component_result, cycle_result, feeder_slot_result, placement_result, head_sequence diff --git a/base_optimizer/smopt_twophase.py b/base_optimizer/smopt_twophase.py index edab124..5ffa528 100644 --- a/base_optimizer/smopt_twophase.py +++ b/base_optimizer/smopt_twophase.py @@ -1,4 +1,6 @@ from base_optimizer.optimizer_common import * +from base_optimizer.smtopt_route import * +from base_optimizer.result_analysis import * def list_range(start, end=None): @@ -80,12 +82,13 @@ def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, part objbst, objbnd = mdl.cbGet(GRB.Callback.MIP_OBJBST), mdl.cbGet(GRB.Callback.MIP_OBJBND) changetime = mdl.cbGet(GRB.Callback.RUNTIME) nonlocal pre_objbst, pre_changetime - # condition: value change - if pre_objbst and abs(pre_objbst - objbst) < 1e-3: - if pre_changetime and changetime - pre_changetime > 100 * (1 - objbnd / objbst): - mdl.terminate() - else: - pre_changetime = changetime +# condition: value change + if abs(objbst - 1e+100) > 1: # 避免未找到可行解提前退出 + if pre_objbst and abs(pre_objbst - objbst) < 1e-3: + if pre_changetime and changetime - pre_changetime > 90 * (1 - objbnd / objbst): + mdl.terminate() + else: + pre_changetime = changetime pre_objbst = objbst @@ -147,7 +150,7 @@ def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, part level += 1 L = len(cycle_assignment) if partition else len(pcb_data) - S = ratio * sum(component_feeder.values()) * 2 if len(feeder_data) == 0 else arg_slot_rng[-1] - arg_slot_rng[0] + 1 # the available feeder num + S = ratio * sum(component_feeder.values()) if len(feeder_data) == 0 else arg_slot_rng[-1] - arg_slot_rng[0] + 1 # the available feeder num M = len(pcb_data) # a sufficiently large number (number of placement points) HC = [[0 for _ in range(J)] for _ in range(I)] for i in range(I): @@ -157,12 +160,12 @@ def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, part mdl = Model('SMT') mdl.setParam('Seed', 0) mdl.setParam('OutputFlag', hinter) # set whether output the debug information - mdl.setParam('TimeLimit', 3600 * 3) + mdl.setParam('TimeLimit', 3600) mdl.setParam('PoolSearchMode', 2) - mdl.setParam('PoolSolutions', 3e2) + mdl.setParam('PoolSolutions', 100) mdl.setParam('PoolGap', 1e-4) # mdl.setParam('MIPFocus', 2) - # mdl.setParam("Heuristics", 0.5) + mdl.setParam("Heuristics", 0.5) # Use only if other methods, including exploring the tree with the default settings, do not yield a viable solution # mdl.setParam("ZeroObjNodes", 100) @@ -178,7 +181,6 @@ def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, part c[i, h, l] <= component_list[cpidx_2_part[i]] for i in range(I) for h in range(max_head_index) for l in range(L)) - # todo: the condition for upper limits of feeders exceed 1 f = {} for i in range(I): if i not in part_feederbase.keys(): @@ -227,7 +229,7 @@ def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, part continue part_feederbase[part_2_cpidx[part]] = slot - # f[slot, part_2_cpidx[part]].Start = 1 + f[slot, part_2_cpidx[part]].Start = 1 slot += 1 for idx, cycle in enumerate(cycle_index): @@ -334,13 +336,12 @@ def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, part in range(L)) if reduction: - # mdl.addConstrs(WL[l] >= WL[l + 1] for l in range(L - 1)) + mdl.addConstrs(WL[l] >= WL[l + 1] for l in range(L - 1)) mdl.addConstr(quicksum(WL[l] for l in range(L)) <= sum(cycle_assignment)) mdl.addConstr(quicksum(WL[l] for l in range(L)) >= math.ceil(len(pcb_data) / max_head_index)) mdl.addConstrs(quicksum(z[j, h, l] for j in range(J) for h in range(max_head_index)) >= quicksum( z[j, h, l + 1] for j in range(J) for h in range(max_head_index)) for l in range(L - 1)) - mdl.addConstrs(y[i, h, l] <= WL[l] for i in range(I) for h in range(max_head_index) for l in range(L)) mdl.addConstrs(v[s, h, l] <= WL[l] for s in range(S) for h in range(max_head_index) for l in range(L)) @@ -365,140 +366,142 @@ def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, part # mdl.optimize() # === result generation === - nozzle_assign, component_assign = [], [] - feeder_assign, cycle_assign = [], [] + opt_res_list = defaultdict(OptResult) if mdl.Status == GRB.OPTIMAL or mdl.Status == GRB.INTERRUPTED or mdl.Status == GRB.TIME_LIMIT: # === selection from solution pool === - component_pos, component_avg_pos = defaultdict(list), defaultdict(list) + component_pos = defaultdict(list[Point]) for _, data in pcb_data.iterrows(): component_index = component_data[component_data.part == data.part].index.tolist()[0] - component_pos[component_index].append([data.x, data.y]) + component_pos[component_index].append(Point(data.x, data.y)) - for i in component_pos.keys(): - component_pos[i] = sorted(component_pos[i], key=lambda pos: (pos[0], pos[1])) - component_avg_pos[i] = [sum(map(lambda pos: pos[0], component_pos[i])) / len(component_pos[i]), - sum(map(lambda pos: pos[1], component_pos[i])) / len(component_pos[i])] + for part in component_pos.keys(): + component_pos[part] = sorted(component_pos[part], key=lambda pos: (pos.x, pos.y)) - min_dist, solution_number = None, -1 + min_dist, solution_number = None, 0 for sol_counter in range(mdl.SolCount): - nozzle_assign, component_assign = [], [] - feeder_assign, cycle_assign = [], [] - mdl.Params.SolutionNumber = sol_counter - pos_counter = defaultdict(int) - - dist = 0 - cycle_placement, cycle_points = defaultdict(list), defaultdict(list) + opt_res = OptResult() + # == 转换标准的贴装头分配的解 === for l in range(L): if abs(WL[l].Xn) <= 1e-4: continue - cycle_placement[l], cycle_points[l] = [-1] * max_head_index, [None] * max_head_index + opt_res.cycle_assign.append(round(WL[l].Xn)) + opt_res.component_assign.append([-1] * max_head_index) + opt_res.feeder_slot_assign.append([-1] * max_head_index) - for h in range(max_head_index): - for l in range(L): - if abs(WL[l].Xn) <= 1e-4: - continue - - pos_list = [] + for h in range(max_head_index): for i in range(I): if abs(y[i, h, l].Xn) <= 1e-4: continue + opt_res.component_assign[-1][h] = i - for _ in range(round(WL[l].Xn)): - pos_list.append(component_pos[i][pos_counter[i]]) - pos_counter[i] += 1 + for s in range(S): + if abs(v[s, h, l].Xn - 1) < 1e-4 and opt_res.component_assign[-1][h] != -1: + opt_res.feeder_slot_assign[-1][h] = s - cycle_placement[l][h] = i - cycle_points[l][h] = [sum(map(lambda pos: pos[0], pos_list)) / len(pos_list), - sum(map(lambda pos: pos[1], pos_list)) / len(pos_list)] - for l in range(L): - if abs(WL[l].Xn) <= 1e-4: - continue + # 根据贴装头位置，转换供料器槽位 + cp_avg_head, cp_sum_cycle = defaultdict(float), defaultdict(int) + for cycle, component_assign in enumerate(opt_res.component_assign): + for head, part in enumerate(component_assign): + if part == -1: + continue + cp_avg_head[part] += opt_res.cycle_assign[cycle] * head + cp_sum_cycle[part] += opt_res.cycle_assign[cycle] - if min_dist is None or dist < min_dist: - min_dist = dist - solution_number = sol_counter + for part, head in cp_avg_head.items(): + cp_avg_head[part] = head / cp_sum_cycle[part] - mdl.Params.SolutionNumber = solution_number - # === 更新吸嘴、元件、周期数优化结果 === - for l in range(L): - nozzle_assign.append([-1 for _ in range(max_head_index)]) - component_assign.append([-1 for _ in range(max_head_index)]) - feeder_assign.append([-1 for _ in range(max_head_index)]) + avg_position = sum([data.x - cp_avg_head[part_2_cpidx[data.part]] * head_interval for _, data in + pcb_data.iterrows()]) / len(pcb_data) + avg_slot = 0 + D_PU, D_PL, D_BW, D_FW = 0, 0, 0, 0 + for cycle, slots in enumerate(opt_res.feeder_slot_assign): + min_slot, max_slot = max_slot_index, 0 + for head, slot in enumerate(slots): + if slot == -1: + continue + min_slot = min(min_slot, slot - head * ratio) + max_slot = max(max_slot, slot - head * ratio) + avg_slot += (max_slot - min_slot) * opt_res.cycle_assign[cycle] + D_PU += (max_slot - min_slot) * slot_interval * opt_res.cycle_assign[cycle] # 拾取路径 - cycle_assign.append(round(WL[l].Xn)) - if abs(WL[l].Xn) <= 1e-4: - continue + avg_slot /= sum(opt_res.cycle_assign) + start_slot = round((avg_position + stopper_pos[0] - slotf1_pos[0]) / slot_interval + avg_slot / 2) + 1 - for h in range(max_head_index): - for i in range(I): - if abs(y[i, h, l].Xn - 1) < 1e-4: - component_assign[-1][h] = i + for cycle in range(len(opt_res.feeder_slot_assign)): + for head in range(max_head_index): + if (slot := opt_res.feeder_slot_assign[cycle][head]) == -1: + continue + opt_res.feeder_slot_assign[cycle][head] = start_slot + slot * (2 if ratio == 1 else 1) - for j in range(J): - if HC[i][j]: - nozzle_assign[-1][h] = j - for s in range(S): - if abs(v[s, h, l].Xn - 1) < 1e-4 and component_assign[l][h] != -1: - feeder_assign[l][h] = s - - # === 更新供料器分配结果 == - component_head = defaultdict(int) - for i in range(I): - cycle_num = 0 - for l, component_cycle in enumerate(component_assign): - for head, component in enumerate(component_cycle): - if component == i: - component_head[i] += cycle_assign[l] * head - cycle_num += cycle_assign[l] - component_head[i] /= cycle_num # 不同元件的加权拾取贴装头 - - average_pos = 0 - for _, data in pcb_data.iterrows(): - average_pos += (data.x - component_head[part_2_cpidx[data.part]] * head_interval) - - average_pos /= len(pcb_data) # 实际贴装位置的加权平均 - average_slot = 0 - for l in range(L): - if abs(WL[l].Xn) <= 1e-4: - continue - min_slot, max_slot = None, None + component_pos_counter = defaultdict(int) + cycle_place_pos = defaultdict(list[Point]) for head in range(max_head_index): - if abs(WL[l].Xn) <= 1e-4 or feeder_assign[l][head] == -1: - continue - slot = feeder_assign[l][head] - head * ratio - if min_slot is None or slot < min_slot: - min_slot = slot - if max_slot is None or slot > max_slot: - max_slot = slot - average_slot += (max_slot - min_slot) * cycle_assign[l] - average_slot /= sum(cycle_assign) - start_slot = round((average_pos + stopper_pos[0] - slotf1_pos[0]) / slot_interval + average_slot / 2) + 1 + for cycle in range(len(opt_res.cycle_assign)): + if (part := opt_res.component_assign[cycle][head]) == -1: + continue + avg_place_pos = Point(0, 0, _h=head) + for counter in range(round(opt_res.cycle_assign[cycle])): + avg_place_pos.x = (1 - 1.0 / (counter + 1)) * avg_place_pos.x + ( + component_pos[part][component_pos_counter[part]].x - head * head_interval) / ( + counter + 1) + avg_place_pos.y = (1 - 1.0 / (counter + 1)) * avg_place_pos.y + component_pos[part][ + component_pos_counter[part]].y / (counter + 1) + component_pos_counter[part] += 1 + avg_place_pos.x += stopper_pos[0] + avg_place_pos.y += stopper_pos[1] + cycle_place_pos[cycle].append(avg_place_pos) + + for cycle in range(len(opt_res.cycle_assign)): + min_slot, max_slot = max_slot_index, 0 + for head in range(max_head_index): + if (slot := opt_res.feeder_slot_assign[cycle][head]) == -1: + continue + min_slot = min(min_slot, slot - head * interval_ratio) + max_slot = max(max_slot, slot - head * interval_ratio) + # cycle_place_pos[cycle] = sorted(cycle_place_pos[cycle], key=lambda pt: pt.x) + + pick_pos = Point(slotf1_pos[0] + (min_slot + max_slot) / 2 * slot_interval, slotf1_pos[1]) + _, seq = dynamic_programming_cycle_route(cycle_place_pos[cycle], pick_pos) + head_position = [Point(0, 0) for _ in range(max_head_index)] + for point in cycle_place_pos[cycle]: + head_position[point.h] = point + + for idx in range(len(seq) - 1): + h1, h2 = seq[idx], seq[idx + 1] + D_PL += max(abs(head_position[h1].x - head_position[h2].x), + abs(head_position[h1].y - head_position[h2].y)) * opt_res.cycle_assign[cycle] + + + # opt_res.placement_assign, opt_res.head_sequence = scan_based_placement_route_generation(component_data, + # pcb_data, + # opt_res.component_assign, + # opt_res.cycle_assign, + # opt_res.feeder_slot_assign, + # hinter=False) + # info = placement_info_evaluation(component_data, pcb_data, opt_res) + # print(f'{info.place_distance + info.pickup_distance: .3f}\t{D_PL + D_PU: .3f}') + opt_res_list[sol_counter] = opt_res + + solution_number = 0 + mdl.Params.SolutionNumber = 0 + if hinter: + print('total cost = {}'.format(mdl.objval)) + print('cycle = {}, nozzle change = {}, pick up = {}'.format(quicksum(WL[l].Xn for l in range(L)), quicksum( + NC[h].Xn for h in range(max_head_index)), quicksum( + PU[s, l].Xn for s in range(-(max_head_index - 1) * ratio, S) for l in range(L)))) + + print('workload: ') for l in range(L): - if abs(WL[l].Xn) <= 1e-4: - continue + print(WL[l].Xn, end=', ') - for h in range(max_head_index): - for s in range(S): - if abs(v[s, h, l].Xn - 1) < 1e-4 and component_assign[l][h] != -1: - feeder_assign[l][h] = start_slot + s * (2 if ratio == 1 else 1) + print('') + print('result') + print('component assignment: ', opt_res_list[solution_number].component_assign) + print('feeder assignment: ', opt_res_list[solution_number].feeder_slot_assign) + print('cycle assignment: ', opt_res_list[solution_number].cycle_assign) - if hinter: - print('total cost = {}'.format(mdl.objval)) - print('cycle = {}, nozzle change = {}, pick up = {}'.format(quicksum(WL[l].Xn for l in range(L)), quicksum( - NC[h].Xn for h in range(max_head_index)), quicksum( - PU[s, l].Xn for s in range(-(max_head_index - 1) * ratio, S) for l in range(L)))) - - print('workload: ') - for l in range(L): - print(WL[l].Xn, end=', ') - - print('') - print('result') - print('nozzle assignment: ', nozzle_assign) - print('component assignment: ', component_assign) - print('feeder assignment: ', feeder_assign) - print('cycle assignment: ', cycle_assign) - return component_assign, feeder_assign, cycle_assign + return opt_res_list[solution_number].component_assign, opt_res_list[solution_number].feeder_slot_assign, \ + opt_res_list[solution_number].cycle_assign diff --git a/base_optimizer/smtopt_route.py b/base_optimizer/smtopt_route.py index 279bad5..6c793b1 100644 --- a/base_optimizer/smtopt_route.py +++ b/base_optimizer/smtopt_route.py @@ -1,3 +1,12 @@ +import copy +import math +import random +import time + +# import matplotlib.pyplot as plt +import numpy as np +import openpyxl + from base_optimizer.optimizer_common import * @@ -13,10 +22,6 @@ def dynamic_programming_cycle_path(pcb_data, cycle_placement, assigned_feeder): head_set.append(head) placement = cycle_placement[head] - if feeder != -1 and placement == -1: - print(assigned_feeder) - print(cycle_placement) - pos.append([pcb_data.iloc[placement]['x'] - head * head_interval + stopper_pos[0], pcb_data.iloc[placement]['y'] + stopper_pos[1], pcb_data.iloc[placement]['r'], head]) @@ -82,9 +87,101 @@ def dynamic_programming_cycle_path(pcb_data, cycle_placement, assigned_feeder): return ans_dist, head_sequence +def dynamic_programming_cycle_route(place_pos: list[Point], pick_pos: Point, rtn_seq=True): + head_sequence = [] + num_pos = len(place_pos) + 1 + + pos = [pick_pos] + place_pos + + def get_distance(pos_1: Point, pos_2: Point): + # 拾取起始与终止位置 或 非同轴 + if pos_1.h is None or pos_2.h is None or pos_1.h + (1 if pos_1.h % 2 == 0 else -1) != pos_2.h: + return max(axis_moving_time(pos_1.x - pos_2.x, 0), axis_moving_time(pos_1.y - pos_2.y, 1)) + else: + return max(axis_moving_time(pos_1.x - pos_2.x, 0), axis_moving_time(pos_1.y - pos_2.y, 1), + head_rotary_time(pos_1.r - pos_2.r)) + + # 各节点之间的距离 + dist = [[get_distance(pos_1, pos_2) for pos_2 in pos] for pos_1 in pos] + + min_dist = [[np.inf for _ in range(num_pos)] for s in range(1 << num_pos)] + min_path = [[[] for _ in range(num_pos)] for s in range(1 << num_pos)] if rtn_seq else None + + # 状压dp搜索 + for s in range(1, 1 << num_pos, 2): + # 考虑节点集合s必须包括节点0 + if not (s & 1): + continue + for j in range(1, num_pos): + # 终点j需在当前考虑节点集合s内 + if not (s & (1 << j)): + continue + if s == int((1 << j) | 1): + # 若考虑节点集合s仅含节点0和节点j，dp边界，赋予初值 + # print('j:', j) + if min_path: + min_path[s][j] = [j] + min_dist[s][j] = dist[0][j] + + # 枚举下一个节点i，更新 + for i in range(1, num_pos): + # 下一个节点i需在考虑节点集合s外 + if s & (1 << i): + continue + if min_dist[s][j] + dist[j][i] < min_dist[s | (1 << i)][i]: + if min_path: + min_path[s | (1 << i)][i] = min_path[s][j] + [i] + min_dist[s | (1 << i)][i] = min_dist[s][j] + dist[j][i] + + ans_dist = float('inf') + ans_path = [] + # 求最终最短哈密顿回路 + for i in range(1, num_pos): + if min_dist[(1 << num_pos) - 1][i] + dist[i][0] < ans_dist: + # 更新，回路化 + if min_path: + ans_path = min_path[s][i] + ans_dist = min_dist[(1 << num_pos) - 1][i] + dist[i][0] + + if len(ans_path): + for parent in ans_path: + head_sequence.append(place_pos[parent - 1].h) + + if place_pos[ans_path[0] - 1].x > place_pos[ans_path[-1] - 1].x: + head_sequence = list(reversed(head_sequence)) + + return ans_dist, head_sequence + + +def quick_sort_cycle_route(place_pos: list[Point], pick_pos: Point, rtn_seq=True): + cycle_place_pos = sorted(place_pos, key=lambda pt: pt.x) + + move_time = 0 + for mount_seq in range(len(cycle_place_pos) - 1): + if cycle_place_pos[mount_seq].h + (1 if cycle_place_pos[mount_seq].h % 2 == 0 else -1) \ + != cycle_place_pos[mount_seq + 1].h: + move_time += max(axis_moving_time(cycle_place_pos[mount_seq].x - cycle_place_pos[mount_seq + 1].x, 0), + axis_moving_time(cycle_place_pos[mount_seq].y - cycle_place_pos[mount_seq + 1].y, 1)) + else: + move_time += max(axis_moving_time(cycle_place_pos[mount_seq].x - cycle_place_pos[mount_seq + 1].x, 0), + axis_moving_time(cycle_place_pos[mount_seq].y - cycle_place_pos[mount_seq + 1].y, 1), + head_rotary_time(cycle_place_pos[mount_seq].r - cycle_place_pos[mount_seq + 1].r)) + + delta_x = axis_moving_time(cycle_place_pos[0].x - pick_pos.x, 0) + delta_y = axis_moving_time(cycle_place_pos[0].y - pick_pos.y, 1) + move_time += max(delta_x, delta_y) + 0.01 * (delta_x ** 2 + delta_y ** 2) + + delta_x = axis_moving_time(cycle_place_pos[-1].x - pick_pos.x, 0) + delta_y = axis_moving_time(cycle_place_pos[-1].y - pick_pos.y, 1) + move_time += max(delta_x, delta_y) + 0.01 * (delta_x ** 2 + delta_y ** 2) + + head_seq = [place_pos.h for place_pos in cycle_place_pos] + return move_time, head_seq + + @timer_wrapper -def greedy_placement_route_generation(component_data, pcb_data, component_result, cycle_result, feeder_slot_result, - hinter=True): +def place_allocate_sequence_route_generation(component_data, pcb_data, component_result, cycle_result, + feeder_slot_result, hinter=True): placement_result, head_sequence_result = [], [] if len(pcb_data) == 0: return placement_result, head_sequence_result @@ -98,7 +195,7 @@ def greedy_placement_route_generation(component_data, pcb_data, component_result mount_point_index[component_index].append(i) mount_point_pos[component_index].append([pcb_data.iloc[i]['x'], pcb_data.iloc[i]['y']]) - search_dir = 1 # 0：自左向右搜索 1：自右向左搜索 + search_dir = 0 # 0：自左向右搜索 1：自右向左搜索 for cycle_set in range(len(component_result)): floor_cycle, ceil_cycle = sum(cycle_result[:cycle_set]), sum(cycle_result[:(cycle_set + 1)]) for cycle in range(floor_cycle, ceil_cycle): @@ -118,7 +215,7 @@ def greedy_placement_route_generation(component_data, pcb_data, component_result continue component_index = component_result[cycle_set][head] - if way_point is None or head_counter % point2head_range == 0: + if way_point is None: index = 0 if way_point is None: if search_dir: @@ -134,6 +231,7 @@ def greedy_placement_route_generation(component_data, pcb_data, component_result [abs(mount_point_pos[component_index][i][0] - way_point[0]) * .1 + abs( mount_point_pos[component_index][i][1] - way_point[1]) for i in range(num_points)]) + # index = random.randint(0, num_points - 1) head = next_head break # index = np.argmax(mount_point_pos[component_index], axis=0)[0] @@ -173,7 +271,7 @@ def greedy_placement_route_generation(component_data, pcb_data, component_result component_index = component_result[cycle_set][head_index] assert 0 <= head_index < max_head_index - + # point_index = random.randint(0, len(mount_point_index[component_index]) - 1) # 贴装点随机分配 assigned_placement[head_index] = mount_point_index[component_index][point_index] way_point = mount_point_pos[component_index][point_index] way_point[0] += (max_head_index - head_index - 1) * head_interval if search_dir \ @@ -191,7 +289,7 @@ def greedy_placement_route_generation(component_data, pcb_data, component_result @timer_wrapper def beam_search_route_generation(component_data, pcb_data, component_result, cycle_result, feeder_slot_result): - beam_width = 10 # 集束宽度 + beam_width = 4 # 集束宽度 base_points = [float('inf'), float('inf')] mount_point_index = [[] for _ in range(len(component_data))] @@ -205,10 +303,8 @@ def beam_search_route_generation(component_data, pcb_data, component_result, cyc mount_point_pos[component_index].append([data.x, data.y]) # 记录最左下角坐标 - if mount_point_pos[component_index][-1][0] < base_points[0]: - base_points[0] = mount_point_pos[component_index][-1][0] - if mount_point_pos[component_index][-1][1] < base_points[1]: - base_points[1] = mount_point_pos[component_index][-1][1] + base_points[0] = min(base_points[0], mount_point_pos[component_index][-1][0]) + base_points[1] = min(base_points[1], mount_point_pos[component_index][-1][1]) beam_placement_sequence, beam_head_sequence = [], [] beam_mount_point_index, beam_mount_point_pos = [], [] @@ -365,21 +461,23 @@ def beam_search_route_generation(component_data, pcb_data, component_result, cyc return beam_placement_sequence[index], beam_head_sequence[index] -def scan_based_placement_route_generation(component_data, pcb_data, component_assign, cycle_assign, feeder_slot_result): +# @timer_wrapper +def scan_based_placement_route_generation(component_data, pcb_data, component_assign, cycle_assign, feeder_slot_result, + hinter=True): placement_result, head_sequence_result = [], [] mount_point_pos, mount_point_index, mount_point_angle, mount_point_part = [], [], [], [] - for i, data in pcb_data.iterrows(): + for _, data in pcb_data.iterrows(): component_index = component_data[component_data.part == data.part].index.tolist()[0] # 记录贴装点序号索引和对应的位置坐标 - mount_point_index.append(i) + mount_point_index.append(len(mount_point_index)) mount_point_pos.append([data.x + stopper_pos[0], data.y + stopper_pos[1]]) mount_point_angle.append(data.r) mount_point_part.append(component_index) - lBoundary, rBoundary = min(mount_point_pos, key=lambda x: x[0])[0], max(mount_point_pos, key=lambda x: x[0])[0] - search_step = max((rBoundary - lBoundary) / max_head_index / 2, 0) + left_boundary, right_boundary = min(mount_point_pos, key=lambda x: x[0])[0], max(mount_point_pos, key=lambda x: x[0])[0] + search_step = max((right_boundary - left_boundary) / max_head_index / 2, 0) ref_pos_y = min(mount_point_pos, key=lambda x: x[1])[1] for cycle_index, component_cycle in enumerate(component_assign): @@ -390,15 +488,15 @@ def scan_based_placement_route_generation(component_data, pcb_data, component_as for search_dir in range(3): # 不同的搜索方向，贴装头和起始点的选取方法各不相同 if search_dir == 0: # 从左向右搜索 - searchPoints = np.arange(lBoundary, (lBoundary + rBoundary) / 2, search_step) + searchPoints = np.arange(left_boundary, (left_boundary + right_boundary) / 2, search_step) head_range = list(range(max_head_index)) elif search_dir == 1: # 从右向左搜索 - searchPoints = np.arange(rBoundary + 1e-3, (lBoundary + rBoundary) / 2, -search_step) + searchPoints = np.arange(right_boundary + 1e-3, (left_boundary + right_boundary) / 2, -search_step) head_range = list(range(max_head_index - 1, -1, -1)) else: # 从中间向两边搜索 - searchPoints = np.arange(lBoundary, rBoundary, search_step / 2) + searchPoints = np.arange(left_boundary, right_boundary, search_step / 2) head_range, head_index = [], (max_head_index - 1) // 2 while head_index >= 0: if 2 * head_index != max_head_index - 1: @@ -426,7 +524,7 @@ def scan_based_placement_route_generation(component_data, pcb_data, component_as for index, mount_index in enumerate(mount_point_index_cpy): if mount_point_part[mount_index] != component_index: continue - horizontal_distance = abs(mount_point_pos_cpy[index][0] - startPoint) + 1e-3 * abs( + horizontal_distance = abs(mount_point_pos_cpy[index][0] - startPoint) + 0 * abs( mount_point_pos_cpy[index][1] - ref_pos_y) if min_horizontal_distance is None or horizontal_distance < min_horizontal_distance: @@ -437,8 +535,8 @@ def scan_based_placement_route_generation(component_data, pcb_data, component_as min_cheby_distance = None next_comp_index = component_assign[cycle_index][head_index] - if assigned_placement[head_index] != -1 or next_comp_index == -1: - continue + # if assigned_placement[head_index] != -1 or next_comp_index == -1: + # continue for index, mount_index in enumerate(mount_point_index_cpy): if mount_point_part[mount_index] != next_comp_index: continue @@ -455,13 +553,12 @@ def scan_based_placement_route_generation(component_data, pcb_data, component_as point_pos = sorted(point_pos, key=lambda x: x[0]) for mount_seq in range(len(point_pos) - 1): - cheby_distance += max(abs(point_pos[mount_seq][0] - point_pos[mount_seq + 1][0]), - abs(point_pos[mount_seq][1] - point_pos[mount_seq + 1][1])) - euler_distance += math.sqrt( - (point_pos[mount_seq][0] - point_pos[mount_seq + 1][0]) ** 2 + ( - point_pos[mount_seq][1] - point_pos[mount_seq + 1][1]) ** 2) + delta_x = axis_moving_time(point_pos[mount_seq][0] - point_pos[mount_seq + 1][0], 0) + delta_y = axis_moving_time(point_pos[mount_seq][1] - point_pos[mount_seq + 1][1], 1) + cheby_distance += max(delta_x, delta_y) + euler_distance += math.sqrt(delta_x ** 2 + delta_y ** 2) - cheby_distance += 0.01 * euler_distance + # cheby_distance += 0.01 * euler_distance if min_cheby_distance is None or cheby_distance < min_cheby_distance: min_cheby_distance, min_euler_distance = cheby_distance, euler_distance point_index = index @@ -687,8 +784,7 @@ def placement_route_relink_heuristic(component_data, pcb_data, placement_result, cycle_average_pos[chg_cycle_index] = [sum(map(lambda x: x[0], point_list)) / len(point_list), sum(map(lambda x: x[1], point_list)) / len(point_list)] - test1 = sum(cycle_length) - test2 = sum(best_cycle_length) + if sum(cycle_length) < sum(best_cycle_length): best_cycle_length = cycle_length.copy() best_cycle_average_pos = copy.deepcopy(cycle_average_pos) @@ -703,4 +799,926 @@ def placement_route_relink_heuristic(component_data, pcb_data, placement_result, prev_time = cur_time # print("number of iteration: ", n_iteration) - return best_placement_result, best_head_sequence_result \ No newline at end of file + return best_placement_result, best_head_sequence_result + + +class RouteNode: + def __init__(self, component_data=None, pcb_data=None): + self.cycle_time = 0 + self.total_time = 0 + self.placement_res = [] + self.headseq_res = [] + + if component_data is not None and pcb_data is not None: + self.unassigned_comp_point = [[] for _ in range(len(component_data))] + for idx, data in pcb_data.iterrows(): + component_index = component_data[component_data['part'] == data.part].index.tolist()[0] + self.unassigned_comp_point[component_index].append(idx) + else: + self.unassigned_comp_point = [[] for _ in range(100)] + + def __lt__(self, other): + return self.cycle_time < other.cycle_time + + +@timer_wrapper +def dynamic_beam_route_generation(component_data, pcb_data, component_assign, cycle_assign, feeder_assign): + empty_node = RouteNode(component_data, pcb_data) + max_beam_width = 4 + point_pos = [Point(data.x + stopper_pos[0], data.y + stopper_pos[1], data.r) for _, data in pcb_data.iterrows()] + + left_boundary, right_boundary = min(point_pos, key=lambda pt: pt.x).x, max(point_pos, key=lambda pt: pt.x).x + search_step = max((right_boundary - left_boundary) / max_head_index, 0) + cycle_pick_pos = [] + for feeder_slot in feeder_assign: + slot_set = set() + for head, slot in enumerate(feeder_slot): + if slot == -1: + continue + slot_set.add(slot - head * interval_ratio) + pick_pos = Point(slotf1_pos[0] + ((min(list(slot_set)) + max(list(slot_set)))/ 2 - 1) * slot_interval, slotf1_pos[1], 0) + cycle_pick_pos.append(pick_pos) + + ref_node = None + placement_result, head_sequence_result = [], [] + start_time = time.time() + for beam_width in range(1, max_beam_width + 1): + beam_node_list = [copy.deepcopy(empty_node) for _ in range(beam_width)] + cycle_index = 0 + current_ref_node = RouteNode(component_data, pcb_data) + for cycle_group_index, cycle_num in enumerate(cycle_assign): + for _ in range(cycle_num): + for beam_node in beam_node_list: + beam_node.placement_res.append([-1 for _ in range(max_head_index)]) + + if ref_node: + current_ref_node.placement_res.append([-1 for _ in range(max_head_index)]) + + beam_node_list_all = [] + for search_dir in range(3): # 不同的搜索方向，贴装头和起始点的选取方法各不相同 + if search_dir == 0: + # 从左向右搜索 + ref_pos_x_list = np.arange(left_boundary, (left_boundary + right_boundary) / 2, search_step) + head_range = list(range(max_head_index)) + elif search_dir == 1: + # 从右向左搜索 + ref_pos_x_list = np.arange(right_boundary + 1e-3, (left_boundary + right_boundary) / 2, -search_step) + head_range = list(range(max_head_index - 1, -1, -1)) + else: + # 从中间向两边搜索 + ref_pos_x_list = np.arange(left_boundary, right_boundary, search_step) + head_range, head_index = [], (max_head_index - 1) // 2 + while head_index >= 0: + if 2 * head_index != max_head_index - 1: + head_range.append(max_head_index - 1 - head_index) + head_range.append(head_index) + head_index -= 1 + + for ref_pos_x in ref_pos_x_list: + beam_node_assigning_list = copy.deepcopy(beam_node_list) + cur_assigning_ref_node = copy.deepcopy(current_ref_node) + is_first_head = True + for head_index in head_range: + if (part_index := component_assign[cycle_group_index][head_index]) == -1: + continue + + if is_first_head: + is_first_head = False + for beam_index, beam_node in enumerate(beam_node_assigning_list): + point_index, min_horizontal_dist = -1, None + for index in beam_node.unassigned_comp_point[part_index]: + horizontal_dist = abs(point_pos[index].x - ref_pos_x) + if min_horizontal_dist is None or horizontal_dist < min_horizontal_dist: + min_horizontal_dist, point_index = horizontal_dist, index + beam_node.placement_res[-1][head_index] = point_index + beam_node.unassigned_comp_point[part_index].remove(point_index) + else: + beam_move_time, beam_indices, beam_point_indices = [], [], [] + for beam_index, beam_node in enumerate(beam_node_assigning_list): + cycle_place_pos = [] + for next_head in range(max_head_index): + if beam_node.placement_res[-1][next_head] == -1: + continue + cycle_place_pos.append(Point(point_pos[beam_node.placement_res[-1][ + next_head]].x - next_head * head_interval, point_pos[ + beam_node.placement_res[-1][next_head]].y, + point_pos[ + beam_node.placement_res[-1][next_head]].r, + next_head)) + cycle_place_pos.append([]) + for point_index in beam_node.unassigned_comp_point[part_index]: + cycle_place_pos[-1] = Point(point_pos[point_index].x - head_index * head_interval, + point_pos[point_index].y, point_pos[point_index].r, + head_index) + + move_time, _ = quick_sort_cycle_route(cycle_place_pos, + cycle_pick_pos[cycle_group_index]) + + beam_move_time.append(move_time) + beam_indices.append(beam_index) + beam_point_indices.append(point_index) + + next_node_assigning_list, assigned_placement_res = [], [] + for index in np.argsort(beam_move_time): + tmp_placement_res = beam_node_assigning_list[beam_indices[index]].placement_res[-1].copy() + tmp_placement_res[head_index] = beam_point_indices[index] + if tmp_placement_res in assigned_placement_res: + continue + + assigned_placement_res.append(tmp_placement_res) + next_node_assigning_list.append( + copy.deepcopy(beam_node_assigning_list[beam_indices[index]])) + next_node_assigning_list[-1].cycle_time = beam_move_time[index] + next_node_assigning_list[-1].unassigned_comp_point[part_index].remove( + beam_point_indices[index]) + next_node_assigning_list[-1].placement_res[-1][head_index] = beam_point_indices[index] + if len(next_node_assigning_list) == beam_width: + break + + beam_node_assigning_list = next_node_assigning_list + + if ref_node: + point_index = ref_node.placement_res[cycle_index][head_index] + cur_assigning_ref_node.placement_res[-1][head_index] = point_index + cur_assigning_ref_node.unassigned_comp_point[part_index].remove(point_index) + beam_node_assigning_list.append(copy.deepcopy(cur_assigning_ref_node)) + + beam_node_list_all.extend(beam_node_assigning_list) + + for beam_node in beam_node_list_all: + beam_node.cycle_time = 0 + cycle_place_pos = [] + for head, index in enumerate(beam_node.placement_res[-1]): + if index == -1: + continue + cycle_place_pos.append( + Point(point_pos[index].x - head * head_interval, point_pos[index].y, point_pos[index].r, + head)) + + cycle_time, headseq = dynamic_programming_cycle_route(cycle_place_pos, + cycle_pick_pos[cycle_group_index]) + beam_node.headseq_res.append(headseq) + beam_node.total_time += cycle_time + + beam_node_list, assigned_placement_res = [], [] + for index in np.argsort([beam_node.total_time for beam_node in beam_node_list_all]): + if beam_node_list_all[index].placement_res in assigned_placement_res: + continue + assigned_placement_res.append(beam_node_list_all[index].placement_res) + beam_node_list.append(beam_node_list_all[index]) + if len(beam_node_list) == beam_width: + break + + if ref_node: + for head_index in range(max_head_index): + if (part_index := component_assign[cycle_group_index][head_index]) == -1: + continue + + point_index = ref_node.placement_res[cycle_index][head_index] + current_ref_node.placement_res[-1][head_index] = point_index + current_ref_node.unassigned_comp_point[part_index].remove(point_index) + + cycle_place_pos = [] + for head, index in enumerate(current_ref_node.placement_res[-1]): + if index == -1: + continue + cycle_place_pos.append( + Point(point_pos[index].x - head * head_interval, point_pos[index].y, point_pos[index].r, + head)) + + cycle_time, _ = dynamic_programming_cycle_route(cycle_place_pos, cycle_pick_pos[cycle_group_index]) + current_ref_node.total_time += cycle_time + current_ref_node.headseq_res.append(ref_node.headseq_res[cycle_index]) + cycle_index += 1 + + min_idx = 0 + for idx in range(1, len(beam_node_list)): + if beam_node_list[min_idx].total_time > beam_node_list[idx].total_time: + min_idx = idx + print( + f"current beam {beam_width}, assembly time: {beam_node_list[min_idx].total_time:.3f}, running time : " + f"{time.time() - start_time:.3f} s") + + ref_node = copy.deepcopy(beam_node_list[min_idx]) + placement_result, head_sequence_result = ref_node.placement_res, ref_node.headseq_res + return placement_result, head_sequence_result + + +class PAPSolution: + def __init__(self, point_pos: list[Point], cycle_res, feeder_slot_res, placement_res, headseq_res): + self.placement_res = placement_res + self.headseq_res = headseq_res + + self.cycle_center_pos = [] + self.cycle_pick_pos = [] + self.cycle_move_time = [] + + self.unassigned_head = [] + self.unassigned_component = [] + + cycle_index = 0 + for cycle_group_index, cycle_group_num in enumerate(cycle_res): + + slot_set = set() + for head, slot in enumerate(feeder_slot_res[cycle_group_index]): + if slot == -1: + continue + slot_set.add(slot - head * interval_ratio) + pick_pos = Point(slotf1_pos[0] + ((min(list(slot_set)) + max(list(slot_set))) / 2 - 1) * slot_interval, + slotf1_pos[1], 0) + for _ in range(cycle_group_num): + self.cycle_pick_pos.append(pick_pos) + + for _ in range(cycle_group_num): + center_pos = Point(0, 0) + cycle_place_pos = [] + for index, head in enumerate(headseq_res[cycle_index]): + cycle_place_pos.append(Point(point_pos[placement_res[cycle_index][head]].x - head * head_interval, + point_pos[placement_res[cycle_index][head]].y)) + cycle_place_pos[-1].h = head + + center_pos.x = (1 - 1 / (index + 1)) * center_pos.x + cycle_place_pos[-1].x / (index + 1) + center_pos.y = (1 - 1 / (index + 1)) * center_pos.y + cycle_place_pos[-1].y / (index + 1) + + self.cycle_center_pos.append(center_pos) + self.cycle_move_time.append( + dynamic_programming_cycle_route(cycle_place_pos, self.cycle_pick_pos[-1], rtn_seq=False)[0]) + self.unassigned_head.append(None) + self.unassigned_component.append(None) + + cycle_index += 1 + + +def all_random_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + selected_cycle = random.sample(range(0, num_of_cycle),num_of_selected_cycle) + + for cycle in selected_cycle: + selected_head = random.sample(destroy_solution.headseq_res[cycle], 1)[0] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def all_worst_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + selected_cycle = np.argpartition(-np.array(destroy_solution.cycle_move_time), kth=num_of_selected_cycle)[ + :num_of_selected_cycle] + for cycle in selected_cycle: + head_derivative = [] + for head in destroy_solution.headseq_res[cycle]: + delta_x = abs(point_pos[destroy_solution.placement_res[cycle][ + head]].x - head * head_interval - destroy_solution.cycle_center_pos[cycle].x) + delta_y = abs( + point_pos[destroy_solution.placement_res[cycle][head]].y - destroy_solution.cycle_center_pos[cycle].y) + head_derivative.append(math.sqrt(delta_x ** 2 + delta_y ** 2 + 1e-8)) + + selected_head = destroy_solution.headseq_res[cycle][np.argmax(head_derivative)] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def all_weighted_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + total_time = sum(solution.cycle_move_time) + p = [cycle_time / total_time for cycle_time in solution.cycle_move_time] + selected_cycle = np.random.choice(range(num_of_cycle), size=num_of_selected_cycle, p=p, replace=False) + + for cycle in selected_cycle: + head_derivative = [] + for head in destroy_solution.headseq_res[cycle]: + delta_x = abs(point_pos[destroy_solution.placement_res[cycle][ + head]].x - head * head_interval - destroy_solution.cycle_center_pos[cycle].x) + delta_y = abs( + point_pos[destroy_solution.placement_res[cycle][head]].y - destroy_solution.cycle_center_pos[cycle].y) + head_derivative.append(math.sqrt(delta_x ** 2 + delta_y ** 2 + 1e-8)) + + selected_head = random.choices(destroy_solution.headseq_res[cycle], weights=head_derivative, k=1)[0] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def weighted_head_random_cycle_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + selected_cycle = random.sample(range(0, num_of_cycle), num_of_selected_cycle) + + for cycle in selected_cycle: + head_derivative = [] + for head in destroy_solution.headseq_res[cycle]: + delta_x = abs(point_pos[destroy_solution.placement_res[cycle][ + head]].x - head * head_interval - destroy_solution.cycle_center_pos[cycle].x) + delta_y = abs( + point_pos[destroy_solution.placement_res[cycle][head]].y - destroy_solution.cycle_center_pos[cycle].y) + head_derivative.append(math.sqrt(delta_x ** 2 + delta_y ** 2 + 1e-8)) + + selected_head = random.choices(destroy_solution.headseq_res[cycle], weights=head_derivative, k=1)[0] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def random_head_worst_cycle_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + selected_cycle = np.argpartition(-np.array(destroy_solution.cycle_move_time), kth=num_of_selected_cycle)[ + :num_of_selected_cycle] + + for cycle in selected_cycle: + selected_head = random.sample(destroy_solution.headseq_res[cycle], 1)[0] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def worst_head_weighted_cycle_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + total_time = sum(solution.cycle_move_time) + p = [cycle_time / total_time for cycle_time in solution.cycle_move_time] + selected_cycle = np.random.choice(range(num_of_cycle), size=num_of_selected_cycle, p=p, replace=False) + + for cycle in selected_cycle: + head_derivative = [] + for head in destroy_solution.headseq_res[cycle]: + delta_x = abs(point_pos[destroy_solution.placement_res[cycle][ + head]].x - head * head_interval - destroy_solution.cycle_center_pos[cycle].x) + delta_y = abs( + point_pos[destroy_solution.placement_res[cycle][head]].y - destroy_solution.cycle_center_pos[cycle].y) + head_derivative.append(math.sqrt(delta_x ** 2 + delta_y ** 2 + 1e-8)) + + selected_head = destroy_solution.headseq_res[cycle][np.argmax(head_derivative)] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def random_head_weighted_cycle_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + total_time = sum(solution.cycle_move_time) + p = [cycle_time / total_time for cycle_time in solution.cycle_move_time] + selected_cycle = np.random.choice(range(num_of_cycle), size=num_of_selected_cycle, p=p, replace=False) + + for cycle in selected_cycle: + selected_head = random.sample(destroy_solution.headseq_res[cycle], 1)[0] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def worst_head_random_cycle_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + selected_cycle = random.sample(range(0, num_of_cycle), num_of_selected_cycle) + + for cycle in selected_cycle: + head_derivative = [] + for head in destroy_solution.headseq_res[cycle]: + delta_x = abs(point_pos[destroy_solution.placement_res[cycle][ + head]].x - head * head_interval - destroy_solution.cycle_center_pos[cycle].x) + delta_y = abs( + point_pos[destroy_solution.placement_res[cycle][head]].y - destroy_solution.cycle_center_pos[cycle].y) + head_derivative.append(math.sqrt(delta_x ** 2 + delta_y ** 2 + 1e-8)) + + selected_head = destroy_solution.headseq_res[cycle][np.argmax(head_derivative)] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def weighted_head_worst_cycle_break(point_pos, comp_of_point, solution, ratio): + destroy_solution = copy.deepcopy(solution) + unassigned_points = [] + num_of_cycle = len(solution.cycle_move_time) + num_of_selected_cycle = round(num_of_cycle * ratio) + selected_cycle = np.argpartition(-np.array(destroy_solution.cycle_move_time), kth=num_of_selected_cycle)[ + :num_of_selected_cycle] + + for cycle in selected_cycle: + head_derivative = [] + for head in destroy_solution.headseq_res[cycle]: + delta_x = abs(point_pos[destroy_solution.placement_res[cycle][ + head]].x - head * head_interval - destroy_solution.cycle_center_pos[cycle].x) + delta_y = abs( + point_pos[destroy_solution.placement_res[cycle][head]].y - destroy_solution.cycle_center_pos[cycle].y) + head_derivative.append(math.sqrt(delta_x ** 2 + delta_y ** 2 + 1e-8)) + + selected_head = random.choices(destroy_solution.headseq_res[cycle], weights=head_derivative, k=1)[0] + + # update center position + head_num = len(destroy_solution.headseq_res[cycle]) - 1 + point_index = destroy_solution.placement_res[cycle][selected_head] + if head_num == 0: + destroy_solution.cycle_center_pos[cycle] = Point(0, 0) + else: + destroy_solution.cycle_center_pos[cycle].x = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].x - (point_pos[point_index].x - selected_head * head_interval) / head_num + destroy_solution.cycle_center_pos[cycle].y = (1 + 1.0 / head_num) * destroy_solution.cycle_center_pos[ + cycle].y - point_pos[point_index].y / head_num + + # destroy operation + selected_point = destroy_solution.placement_res[cycle][selected_head] + destroy_solution.unassigned_head[cycle] = selected_head + destroy_solution.unassigned_component[cycle] = comp_of_point[selected_point] + + destroy_solution.headseq_res[cycle].remove(selected_head) + destroy_solution.placement_res[cycle][selected_head] = -1 + + unassigned_points.append(selected_point) + + return destroy_solution, unassigned_points + + +def random_repair(point_pos, comp_of_point, solution, unassigned_points): + for point in unassigned_points: + unassigned_cycle_list = [] + for cycle, head in enumerate(solution.unassigned_head): + if head is not None and solution.unassigned_component[cycle] == comp_of_point[point]: + unassigned_cycle_list.append(cycle) + + selected_cycle = random.sample(unassigned_cycle_list, 1)[0] + selected_head = solution.unassigned_head[selected_cycle] + + # update center pos、move time and head sequence + cycle_place_pos = [ + Point(point_pos[point].x - selected_head * head_interval, point_pos[point].y, point_pos[point].r, + selected_head)] + for head in solution.headseq_res[selected_cycle]: + cycle_place_pos.append( + Point(point_pos[solution.placement_res[selected_cycle][head]].x - head * head_interval, + point_pos[solution.placement_res[selected_cycle][head]].y, + point_pos[solution.placement_res[selected_cycle][head]].r, head)) + + solution.cycle_move_time[selected_cycle], solution.headseq_res[ + selected_cycle] = dynamic_programming_cycle_route(cycle_place_pos, solution.cycle_pick_pos[selected_cycle]) + + num_head = len(solution.headseq_res[selected_cycle]) + solution.cycle_center_pos[selected_cycle].x = (1 - 1.0 / num_head) * solution.cycle_center_pos[ + selected_cycle].x + (point_pos[point].x - selected_head * head_interval) / num_head + solution.cycle_center_pos[selected_cycle].y = (1 - 1.0 / num_head) * solution.cycle_center_pos[ + selected_cycle].y + point_pos[point].y / num_head + + solution.placement_res[selected_cycle][selected_head] = point + solution.unassigned_head[selected_cycle] = None + solution.unassigned_component[selected_cycle] = None + return solution + + +def greedy_repair(point_pos, comp_of_point, solution, unassigned_points): + for point in unassigned_points: + unassigned_cycle_list, unassigned_cycle_dist = [], [] + for cycle, head in enumerate(solution.unassigned_head): + if head is not None and solution.unassigned_component[cycle] == comp_of_point[point]: + unassigned_cycle_list.append(cycle) + + delta_x = abs(solution.cycle_center_pos[cycle].x - point_pos[point].x + head * head_interval) + delta_y = abs(solution.cycle_center_pos[cycle].y - point_pos[point].y) + + unassigned_cycle_dist.append(math.sqrt(delta_x ** 2 + delta_y ** 2) + 1e-10) + + selected_cycle = unassigned_cycle_list[np.argmin(unassigned_cycle_dist)] + selected_head = solution.unassigned_head[selected_cycle] + + # update center pos、move time and head sequence + cycle_place_pos = [ + Point(point_pos[point].x - selected_head * head_interval, point_pos[point].y, point_pos[point].r, + selected_head)] + for head in solution.headseq_res[selected_cycle]: + cycle_place_pos.append( + Point(point_pos[solution.placement_res[selected_cycle][head]].x - head * head_interval, + point_pos[solution.placement_res[selected_cycle][head]].y, + point_pos[solution.placement_res[selected_cycle][head]].r, head)) + + solution.cycle_move_time[selected_cycle], solution.headseq_res[ + selected_cycle] = dynamic_programming_cycle_route(cycle_place_pos, solution.cycle_pick_pos[selected_cycle]) + + num_head = len(solution.headseq_res[selected_cycle]) + solution.cycle_center_pos[selected_cycle].x = (1 - 1.0 / num_head) * solution.cycle_center_pos[ + selected_cycle].x + (point_pos[point].x - selected_head * head_interval) / num_head + solution.cycle_center_pos[selected_cycle].y = (1 - 1.0 / num_head) * solution.cycle_center_pos[ + selected_cycle].y + point_pos[point].y / num_head + + solution.placement_res[selected_cycle][selected_head] = point + solution.unassigned_head[selected_cycle] = None + solution.unassigned_component[selected_cycle] = None + return solution + + +def weighted_repair(point_pos, comp_of_point, solution, unassigned_points): + for point in unassigned_points: + unassigned_cycle_list, unassigned_cycle_dist = [], [] + for cycle, head in enumerate(solution.unassigned_head): + if head is not None and solution.unassigned_component[cycle] == comp_of_point[point]: + unassigned_cycle_list.append(cycle) + + delta_x = abs(solution.cycle_center_pos[cycle].x - point_pos[point].x + head * head_interval) + delta_y = abs(solution.cycle_center_pos[cycle].y - point_pos[point].y) + + unassigned_cycle_dist.append(math.sqrt(delta_x ** 2 + delta_y ** 2) + 1e-10) + + selected_cycle = random.choices(unassigned_cycle_list, weights=unassigned_cycle_dist, k=1)[0] + selected_head = solution.unassigned_head[selected_cycle] + + # update center pos、move time and head sequence + cycle_place_pos = [ + Point(point_pos[point].x - selected_head * head_interval, point_pos[point].y, point_pos[point].r, + selected_head)] + for head in solution.headseq_res[selected_cycle]: + cycle_place_pos.append( + Point(point_pos[solution.placement_res[selected_cycle][head]].x - head * head_interval, + point_pos[solution.placement_res[selected_cycle][head]].y, + point_pos[solution.placement_res[selected_cycle][head]].r, head)) + + solution.cycle_move_time[selected_cycle], solution.headseq_res[ + selected_cycle] = dynamic_programming_cycle_route(cycle_place_pos, solution.cycle_pick_pos[selected_cycle]) + + num_head = len(solution.headseq_res[selected_cycle]) + solution.cycle_center_pos[selected_cycle].x = (1 - 1.0 / num_head) * solution.cycle_center_pos[ + selected_cycle].x + (point_pos[point].x - selected_head * head_interval) / num_head + solution.cycle_center_pos[selected_cycle].y = (1 - 1.0 / num_head) * solution.cycle_center_pos[ + selected_cycle].y + point_pos[point].y / num_head + + solution.placement_res[selected_cycle][selected_head] = point + solution.unassigned_head[selected_cycle] = None + solution.unassigned_component[selected_cycle] = None + return solution + + +def select_and_apply_destroy_operator(point_pos, comp_of_point, weight, solution, ratio=0.3): + destroy_operator_map = { + 0: all_random_break, + 1: all_worst_break, + 2: all_weighted_break, + 3: weighted_head_random_cycle_break, + 4: random_head_worst_cycle_break, + 5: worst_head_weighted_cycle_break, + 6: random_head_weighted_cycle_break, + 7: worst_head_random_cycle_break, + 8: weighted_head_worst_cycle_break, + } + assert len(destroy_operator_map) == len(weight) + destroy_index, destroy_solution, unassigned_points = -1, solution, [] + roulette = np.array(weight).cumsum() + r = random.uniform(0, max(roulette)) + for index in range(len(destroy_operator_map)): + if roulette[index] >= r: + destroy_index = index + destroy_solution, unassigned_points = destroy_operator_map[index](point_pos, comp_of_point, solution, ratio) + break + return destroy_index, destroy_solution, unassigned_points + + +def select_and_apply_repair_operator(point_pos, comp_of_point, weight, solution, unassigned_points): + repair_operator_map = { + 0: random_repair, + 1: greedy_repair, + 2: weighted_repair, + } + assert len(repair_operator_map) == len(weight) + repair_index, repair_solution = -1, solution + roulette = np.array(weight).cumsum() + r = random.uniform(0, max(roulette)) + for index in range(len(weight)): + if roulette[index] >= r: + repair_index = index + repair_solution = repair_operator_map[index](point_pos, comp_of_point, solution, unassigned_points) + break + + return repair_index, repair_solution + + +@timer_wrapper +def alns_route_reschedule(pcb_data, placement_res, headseq_res, component_res, feeder_slot_res, cycle_res, hinter=True): + point_pos = [Point(data.x + stopper_pos[0], data.y + stopper_pos[1], data.r) for _, data in pcb_data.iterrows()] + comp_of_point = defaultdict(int) + cycle_index = 0 + for cycle_group_index, cycle_group_num in enumerate(cycle_res): + for _ in range(cycle_group_num): + for index, head in enumerate(headseq_res[cycle_index]): + comp_of_point[placement_res[cycle_index][head]] = component_res[cycle_group_index][head] + cycle_index += 1 + + init_temperature, final_temperature = 1, 0.2 + alpha, beta = 0.8, 0.5 + solution = PAPSolution(point_pos, cycle_res, feeder_slot_res, placement_res, headseq_res) + + num_destroy_operators, num_repair_operators = 9, 3 + destroy_weight, repair_weight = [1 for _ in range(num_destroy_operators)], [1 for _ in range(num_repair_operators)] + destroy_score, repair_score = [1 for _ in range(num_destroy_operators)], [1 for _ in range(num_repair_operators)] + destroy_times, repair_times = [0 for _ in range(num_destroy_operators)], [0 for _ in range(num_repair_operators)] + + best_solution = copy.deepcopy(solution) + max_iteration = 2000 + # === 初始化各周期的移动路径长度，中心位置坐标和未分配的贴装头 === + omega_new_global_best, omega_better_than_current, omega_accept, omega_reject = 3, 2, 0.8, 0.2 + best_move_time_list = [] + with tqdm(total=max_iteration) as pbar: + pbar.set_description('adaptive large neighbor search route reschedule') + for _ in range(max_iteration): + temperature = init_temperature + solution = best_solution + while temperature > final_temperature: + destroy_index, destroy_solution, unassigned_points = select_and_apply_destroy_operator(point_pos, + comp_of_point, + destroy_weight, + solution) + + repair_index, repair_solution = select_and_apply_repair_operator(point_pos, comp_of_point, + repair_weight, destroy_solution, + unassigned_points) + if sum(repair_solution.cycle_move_time) <= sum(solution.cycle_move_time): + solution = repair_solution + if sum(repair_solution.cycle_move_time) <= sum(best_solution.cycle_move_time): + best_solution = repair_solution + destroy_score[destroy_index] += omega_new_global_best + repair_score[repair_index] += omega_new_global_best + else: + destroy_score[destroy_index] += omega_better_than_current + repair_score[repair_index] += omega_better_than_current + else: + if random.random() < np.exp( + (sum(solution.cycle_move_time) - sum(repair_solution.cycle_move_time)) / temperature): + solution = repair_solution + destroy_score[destroy_index] += omega_accept + repair_score[repair_index] += omega_accept + else: + destroy_score[destroy_index] += omega_reject + repair_score[repair_index] += omega_reject + + best_move_time_list.append(sum(best_solution.cycle_move_time)) + destroy_times[destroy_index] += 1 + repair_times[repair_index] += 1 + + destroy_weight[destroy_index] = destroy_weight[destroy_index] * beta + (1 - beta) * destroy_score[ + destroy_index] / destroy_times[destroy_index] + + repair_weight[repair_index] = repair_weight[repair_index] * beta + (1 - beta) * repair_score[ + repair_index] / repair_times[repair_index] + + temperature *= alpha + pbar.update(1) + + # plt.plot(range(len(best_move_time_list)), best_move_time_list) + # plt.show() + + print('best move time', best_move_time_list[-1]) + workbook = openpyxl.load_workbook('result/alns_route.xlsx') + worksheet = workbook.active + writing_colum = None + for col in worksheet.iter_cols(min_col=0, max_col=100): + for cell in col: + if cell.value is None: + writing_colum = cell.column_letter + break + if writing_colum: + break + + for row_num, value in enumerate(best_move_time_list): + cell_reference = f"{writing_colum}{row_num + 1}" + worksheet[cell_reference] = value + workbook.save('result/alns_route.xlsx') + + return best_solution.placement_res, best_solution.headseq_res + + +def place_cluster_greedy_route_generation(component_data, pcb_data, component_result, cycle_result, feeder_slot_result): + placement_result, head_sequence_result = [], [] + + # === assign CT group to feeder slot === + component_point_pos = defaultdict(list) + + for idx, data in pcb_data.iterrows(): + component_point_pos[data.part].append([data.x + stopper_pos[0], data.y + stopper_pos[1], idx]) + + for pos_list in component_point_pos.values(): + pos_list.sort(key=lambda x: (x[0], x[1])) + + component_point_index = defaultdict(int) + for cycle_set in range(len(cycle_result)): + for cycle in range(cycle_result[cycle_set]): + placement_result.append([-1 for _ in range(max_head_index)]) + for head in range(max_head_index): + part_index = component_result[cycle_set][head] + if part_index == -1: + continue + + part = component_data.iloc[part_index]['part'] + point_info = component_point_pos[part][component_point_index[part]] + + placement_result[-1][head] = point_info[2] + # mount_point[head] = point_info[0:2] + + component_point_index[part] += 1 + head_sequence_result.append( + dynamic_programming_cycle_path(pcb_data, placement_result[-1], feeder_slot_result[cycle_set])[1]) + return placement_result, head_sequence_result + + +def greedy_level_placing_route_generation(component_data, pcb_data, component_result, cycle_result, feeder_slot_result, hinter=False): + placement_result, head_sequence = [], [] + part_indices = defaultdict(int) + for part_idx, data in component_data.iterrows(): + part_indices[data.part] = part_idx + + mount_point_pos = defaultdict(list) + for pcb_idx, data in pcb_data.iterrows(): + mount_point_pos[part_indices[data.part]].append([data.x, data.y, pcb_idx]) + + for index_ in mount_point_pos.keys(): + mount_point_pos[index_].sort(key=lambda x: (x[1], x[0])) + + for cycle_idx, _ in enumerate(cycle_result): + for _ in range(cycle_result[cycle_idx]): + placement_result.append([-1 for _ in range(max_head_index)]) + for head in range(max_head_index): + if component_result[cycle_idx][head] == -1: + continue + index_ = component_result[cycle_idx][head] + placement_result[-1][head] = mount_point_pos[index_][-1][2] + mount_point_pos[index_].pop() + head_sequence.append( + dynamic_programming_cycle_path(pcb_data, placement_result[-1], feeder_slot_result[cycle_idx])[1]) + return placement_result, head_sequence + diff --git a/dataloader.py b/dataloader.py index c96f387..d1cdbe9 100644 --- a/dataloader.py +++ b/dataloader.py @@ -161,8 +161,27 @@ def load_data(filename: str, load_feeder=False, auto_register=True): feeder_data = pd.DataFrame(columns=feeder_columns) else: feeder_data = pd.DataFrame(columns=feeder_columns) + + for idx, data in feeder_data.iterrows(): + feeder_data.at[idx, 'slot'] = int(data['slot'][1:]) + feeder_data.sort_values(by='slot', ascending=True, inplace=True, ignore_index=True) + + feeder_data_check = defaultdict(str) + for _, data in feeder_data.iterrows(): + feeder_data_check[data.slot] = data.part + + for machine_index in range(machine_num): + for idx, data in pcb_data[machine_index].iterrows(): + slot, part = data.fdr.split(' ') + slot = int(slot[1:]) + if feeder_data_check[slot] == '': + feeder_data_check[slot] = part + if feeder_data_check[slot] != part: + warning_info = f'conflict feeder registration PCB: {data.fdr}, BASE: F{slot} {feeder_data_check[slot]}' + # warnings.warn(warning_info, UserWarning) + return pcb_data, component_data, feeder_data diff --git a/estimator.py b/estimator.py index 73c5bff..9128238 100644 --- a/estimator.py +++ b/estimator.py @@ -6,9 +6,9 @@ def exact_assembly_time(pcb_data, component_data): feeder_data = pd.DataFrame(columns=['slot', 'part']) component_result, cycle_result, feeder_slot_result = feeder_priority_assignment(component_data, pcb_data, feeder_data, hinter=False) - placement_result, head_sequence_result = greedy_placement_route_generation(component_data, pcb_data, - component_result, cycle_result, - feeder_slot_result, hinter=False) + placement_result, head_sequence_result = place_allocate_sequence_route_generation(component_data, pcb_data, + component_result, cycle_result, + feeder_slot_result, hinter=False) opt_res = OptResult(component_result, cycle_result, feeder_slot_result, placement_result, head_sequence_result) info = placement_info_evaluation(component_data, pcb_data, opt_res) # return info.metric() @@ -19,7 +19,7 @@ def error_info(pred_val, real_val, type='train'): absolute_error = np.array([]) for idx, (t1, t2) in enumerate(np.nditer([pred_val, real_val])): absolute_error = np.append(absolute_error, abs(t1 - t2) / (t2 + 1e-10) * 100) - if absolute_error[-1] > 15: + if absolute_error[-1] > 10: print(f'\033[0;31;31midx: {idx + 1: d}, net: {t1: .3f}, real: {t2: .3f}, ' f'gap: {absolute_error[-1]: .3f}\033[0m') @@ -561,6 +561,8 @@ class MetricEstimator(Estimator): def training(self, params): x_fit, y_fit = data_mgr.metric('opt/' + params.train_file) self.lr.fit(x_fit, y_fit) + y_predict = self.lr.predict(x_fit) + error_info(y_fit, y_predict, 'train') print(self.lr.coef_) def testing(self, params): @@ -578,7 +580,7 @@ if __name__ == '__main__': parser = argparse.ArgumentParser(description='network training implementation') # parser.add_argument('--train', default=True, type=bool, help='determine whether training the network') - parser.add_argument('--save', default=True, type=bool, + parser.add_argument('--save', default=False, type=bool, help='determine whether saving the parameters of network, linear regression model, etc.') parser.add_argument('--overwrite', default=False, type=bool, help='determine whether overwriting the training and testing data') @@ -588,7 +590,7 @@ if __name__ == '__main__': parser.add_argument('--batch_size', default=2000, type=int, help='size of training batch') parser.add_argument('--lr', default=1e-5, type=float, help='learning rate for the network') parser.add_argument('--model', default='neural-network', help='method for assembly time estimation') - parser.add_argument('--machine_optimizer', default='feeder-scan', type=str, help='optimizer for single machine') + parser.add_argument('--machine_optimizer', default='feeder-priority', type=str, help='optimizer for single machine') params = parser.parse_args() data_mgr = DataMgr() @@ -606,13 +608,13 @@ if __name__ == '__main__': # data_mgr.remover() # remove the last saved data # data_mgr.saver('data/' + file_name, pcb_data) # save new data - info = base_optimizer(1, pcb_data, component_data, pd.DataFrame(columns=['slot', 'part', 'arg']), - params, hinter=True) + info = base_optimizer(1, pcb_data, component_data, pd.DataFrame(columns=['slot', 'part']), params, + hinter=True) data_mgr.recorder(f, info, pcb_data, component_data) f.close() - estimator = MetricEstimator() + estimator = NeuralEstimator() estimator.training(params) estimator.testing(params) diff --git a/generator.py b/generator.py index b341789..6708e6b 100644 --- a/generator.py +++ b/generator.py @@ -73,9 +73,9 @@ class DataMgr: lineinfo += '\t' + '{:.3f}'.format(pcb_data['x'].max() - pcb_data['x'].min()) + '\t' + '{:.3f}'.format( pcb_data['y'].max() - pcb_data['y'].min()) - # part_position = defaultdict(list) - # for _, data in pcb_data.iterrows(): - # part_position[data['part']].append([data['x'], data['y']]) + part_position = defaultdict(list) + for _, data in pcb_data.iterrows(): + part_position[data['part']].append([data['x'], data['y']]) point_counter, component_counter = 0, 0 nozzle_type = set() @@ -94,8 +94,8 @@ class DataMgr: if data.points == 0: continue lineinfo += '\t' + data.part + '\t' + data.nz + '\t' + str(data.points) - # lineinfo += '\t' + '{:.3f}'.format(np.ptp([pos[0] for pos in part_position[data.part]])) - # lineinfo += '\t' + '{:.3f}'.format(np.ptp([pos[1] for pos in part_position[data.part]])) + lineinfo += '\t' + '{:.3f}'.format(np.ptp([pos[0] for pos in part_position[data.part]])) + lineinfo += '\t' + '{:.3f}'.format(np.ptp([pos[1] for pos in part_position[data.part]])) lineinfo += '\n' file_handle.write(lineinfo) @@ -173,7 +173,7 @@ class DataMgr: return data def heuristic_objective(self, cp_points, cp_nozzle): - if len(cp_points.keys()) == 0: + if len(cp_points.keys()) or sum(cp_points.values()) == 0: return 0, 0, 0, 0 nozzle_heads, nozzle_points = defaultdict(int), defaultdict(int) for idx, points in cp_points.items(): @@ -417,18 +417,18 @@ class DataMgr: # assembly time data time_data.append(float(items[0])) - cp_data.append([cp_points, cp_nozzle, board_width, board_height]) + cp_data.append([cp_points, cp_nozzle, board_width, board_height, cycle, nozzle_change_data, pickup_data]) line = file.readline() if data_filter: - cph_data = [point_data[idx] / time_data[idx] * 3600 for idx in range(len(time_data))] + # cph_data = [point_data[idx] / time_data[idx] * 3600 for idx in range(len(time_data))] - w_quart = 0.6 - Q1, Q3 = np.percentile(np.array(cph_data), 25), np.percentile(np.array(cph_data), 75) + w_quart = 0.3 + Q1, Q3 = np.percentile(np.array(time_data), 25), np.percentile(np.array(time_data), 75) - indices = [i for i in range(len(cph_data)) if - Q1 - w_quart * (Q3 - Q1) <= cph_data[i] <= Q3 + w_quart * (Q3 - Q1)] + indices = [i for i in range(len(time_data)) if + Q1 - w_quart * (Q3 - Q1) <= time_data[i] <= Q3 + w_quart * (Q3 - Q1)] filter_cp_data, filter_time_data = [], [] for idx in indices: @@ -465,8 +465,8 @@ class DataMgr: def neural_encode(self, input_data): train_data = [] - for cp_points, cp_nozzle, board_width, board_height in input_data: - train_data.append(self.encode(cp_points, cp_nozzle, board_width, board_height)) + for cp_points, cp_nozzle, board_width, board_height, cy, nz, pu in input_data: + train_data.append(self.encode(cp_points, cp_nozzle, board_width, board_height, cy, nz, pu)) return train_data def get_feature(self): diff --git a/lineopt_hyperheuristic.py b/lineopt_hyperheuristic.py index 5ee74c1..0b9eb9a 100644 --- a/lineopt_hyperheuristic.py +++ b/lineopt_hyperheuristic.py @@ -278,7 +278,7 @@ def cal_individual_val(heuristic_map, cp_index, cp_points, cp_nozzle, cp_feeders def line_optimizer_hyperheuristic(component_data, pcb_data, machine_number): heuristic_map = { 'p': LeastPoints, - 'n': LeastNzChange, + 'n': LeastNzTypes, 'c': LeastCpTypes, 'r': LeastCpNzRatio, 'k': LeastCycle, @@ -411,12 +411,13 @@ def line_optimizer_hyperheuristic(component_data, pcb_data, machine_number): best_component_list = component_list.copy() machine_cp_points = convert_assignment_result(heuristic_map, cp_index, cp_points, cp_nozzle, cp_feeders, - best_component_list, best_heuristic_list, machine_number) + best_component_list, best_heuristic_list, machine_number, is_opt=True) assignment_result = [[0 for _ in range(len(component_data))] for _ in range(machine_number)] for machine_idx in range(machine_number): for idx in machine_cp_points[machine_idx]: assignment_result[machine_idx][cp_index[idx]] += cp_points[idx] + return assignment_result diff --git a/lineopt_model.py b/lineopt_model.py index 4472466..5e29f7d 100644 --- a/lineopt_model.py +++ b/lineopt_model.py @@ -7,7 +7,7 @@ def line_optimizer_model(component_data, pcb_data, machine_num, hinter=True): mdl = Model('pcb assembly line optimizer') mdl.setParam('Seed', 0) mdl.setParam('OutputFlag', hinter) # set whether output the debug information - mdl.setParam('TimeLimit', 600 * 3) + mdl.setParam('TimeLimit', 1000) nozzle_type, component_type = [], [] for _, data in component_data.iterrows(): @@ -22,9 +22,10 @@ def line_optimizer_model(component_data, pcb_data, machine_num, hinter=True): H = max_head_index I = len(component_data) - S = min(len(component_data) * ratio, 60) - K = math.ceil(len(pcb_data) * 1.0 / H / M) + 1 - # K = 3 + S = sum([data.fdn * ratio for _, data in component_data.iterrows()]) + K = math.ceil(len(pcb_data) * 1.0 / M) + 1 + + # K = len(pcb_data) CompOfNozzle = [[0 for _ in range(J)] for _ in range(I)] # Compatibility component_point = [0 for _ in range(I)] @@ -35,15 +36,17 @@ def line_optimizer_model(component_data, pcb_data, machine_num, hinter=True): # objective related g = mdl.addVars(list_range(K), list_range(M), vtype=GRB.BINARY) - d = mdl.addVars(list_range(K), list_range(H), list_range(M), lb=0, vtype=GRB.CONTINUOUS) + d = mdl.addVars(list_range(K), list_range(H), list_range(M), vtype=GRB.INTEGER) u = mdl.addVars(list_range(I), list_range(K), list_range(H), list_range(M), vtype=GRB.BINARY) v = mdl.addVars(list_range(S), list_range(K), list_range(H), list_range(M), vtype=GRB.BINARY) - d_plus = mdl.addVars(list_range(J), list_range(K), list_range(H), list_range(M), lb=0, vtype=GRB.CONTINUOUS) - d_minus = mdl.addVars(list_range(J), list_range(K), list_range(H), list_range(M), lb=0, vtype=GRB.CONTINUOUS) - w = mdl.addVars(list_range(K), list_range(M), vtype=GRB.CONTINUOUS) + w = mdl.addVars(list_range(J), list_range(K), list_range(H), list_range(M), vtype=GRB.BINARY) + d_plus = mdl.addVars(list_range(J), list_range(K), list_range(H), list_range(M), vtype=GRB.CONTINUOUS) + d_minus = mdl.addVars(list_range(J), list_range(K), list_range(H), list_range(M), vtype=GRB.CONTINUOUS) + z = mdl.addVars(list_range(K), list_range(M), vtype=GRB.INTEGER) e = mdl.addVars(list_range(-(H - 1) * ratio, S), list_range(K), list_range(M), vtype=GRB.BINARY) f = mdl.addVars(list_range(S), list_range(I), list_range(M), vtype=GRB.BINARY, name='') + t = mdl.addVars(list_range(M), lb=0, ub=N, vtype=GRB.CONTINUOUS) obj = mdl.addVar(lb=0, ub=N, vtype=GRB.CONTINUOUS) mdl.addConstrs(g[k, m] >= g[k + 1, m] for k in range(K - 1) for m in range(M)) @@ -52,24 +55,26 @@ def line_optimizer_model(component_data, pcb_data, machine_num, hinter=True): quicksum(u[i, k, h, m] for i in range(I)) <= g[k, m] for k in range(K) for h in range(H) for m in range(M)) # nozzle no more than 1 for head h and cycle k - mdl.addConstrs(quicksum(CompOfNozzle[i][j] * u[i, k, h, m] for i in range(I) for j in range(J)) <= 1 for k - in range(K) for h in range(H) for m in range(M)) + mdl.addConstrs(quicksum(w[j, k, h, m] for j in range(J)) <= 1 for k in range(K) for h in range(H) for m in range(M)) # work completion mdl.addConstrs( quicksum(u[i, k, h, m] for k in range(K) for h in range(H) for m in range(M)) == component_point[i] for i in range(I)) + # nozzle change - mdl.addConstrs(quicksum(CompOfNozzle[i][j] * u[i, k, h, m] for i in range(I)) - quicksum( - CompOfNozzle[i][j] * u[i, k + 1, h, m] for i in range(I)) == d_plus[j, k, h, m] - d_minus[j, k, h, m] for k in + mdl.addConstrs( + u[i, k, h, m] <= quicksum(CompOfNozzle[i][j] * w[j, k, h, m] for j in range(J)) for i in range(I) for k in + range(K) for h in range(H) for m in range(M)) + + mdl.addConstrs(w[j, k, h, m] - w[j, k + 1, h, m] == d_plus[j, k, h, m] - d_minus[j, k, h, m] for k in range(K - 1) for j in range(J) for h in range(H) for m in range(M)) - mdl.addConstrs(quicksum(CompOfNozzle[i][j] * u[i, K - 1, h, m] for i in range(I)) - quicksum( - CompOfNozzle[i][j] * u[i, 0, h, m] for i in range(I)) == d_plus[j, K - 1, h, m] - d_minus[j, K - 1, h, m] for j - in range(J) for h in range(H) for m in range(M)) + mdl.addConstrs(w[j, 0, h, m] - w[j, K - 1, h, m] == d_plus[j, K - 1, h, m] - d_minus[j, K - 1, h, m] + for j in range(J) for h in range(H) for m in range(M)) - mdl.addConstrs(2 * d[k, h, m] == quicksum(d_plus[j, k, h, m] for j in range(J)) + quicksum( - d_minus[j, k, h, m] for j in range(J)) - 1 for k in range(K) for h in range(H) for m in range(M)) + mdl.addConstrs(d[k, h, m] == quicksum(d_plus[j, k, h, m] for j in range(J)) + quicksum( + d_minus[j, k, h, m] for j in range(J)) for k in range(K) for h in range(H) for m in range(M)) # simultaneous pick for s in range(-(H - 1) * ratio, S): @@ -99,25 +104,29 @@ def line_optimizer_model(component_data, pcb_data, machine_num, hinter=True): for s in range(S) for m in range(M)) # pickup movement - mdl.addConstrs(w[k, m] >= s1 * e[s1, k, m] - s2 * e[s2, k, m] + N * (e[s1, k, m] + e[s2, k, m] - 2) for s1 in + mdl.addConstrs(z[k, m] >= s1 * e[s1, k, m] - s2 * e[s2, k, m] + N * (e[s1, k, m] + e[s2, k, m] - 2) for s1 in range(-(H - 1) * ratio, S) for s2 in range(-(H - 1) * ratio, S) for k in range(K) for m in range(M)) # objective - mdl.addConstrs(obj >= Fit_cy * quicksum(g[k, m] for k in range(K)) + Fit_nz * 2 * quicksum( - d[k, h, m] for h in range(H) for k in range(K)) + Fit_pu * quicksum( - e[s, k, m] for s in range(-(H - 1) * ratio, S) for k in range(K)) + Fit_pl * quicksum( - u[i, k, h, m] for i in range(I) for k in range(K) for h in range(H)) + Fit_mv * head_interval * quicksum( - w[k, m] for k in range(K)) for m in range(M)) + mdl.addConstrs(t[m] == Fit_cy * quicksum(g[k, m] for k in range(K)) + Fit_nz * quicksum( + d[k, h, m] for h in range(H) for k in range(K)) + Fit_pu * quicksum( + e[s, k, m] for s in range(-(H - 1) * ratio, S) for k in range(K)) + Fit_pl * quicksum( + u[i, k, h, m] for i in range(I) for k in range(K) for h in range(H)) + Fit_mv * quicksum( + z[k, m] for k in range(K)) for m in range(M)) + for m in range(M - 1): + mdl.addConstr(t[m] >= t[m + 1]) + + mdl.addConstrs(obj >= t[m] for m in range(M)) mdl.setObjective(obj, GRB.MINIMIZE) - mdl.optimize() + for m in range(M): print(f'machine {m} : cycle : {sum(g[k, m].x for k in range(K))}, ' f'nozzle change : {sum(d[k, h, m].x for h in range(H) for k in range(K))}, ' f'pick up : {sum(e[s, k, m].x for s in range(-(H - 1) * ratio, S) for k in range(K))}, ' f'placement : {sum(u[i, k, h, m].x for i in range(I) for k in range(K) for h in range(H))}, ' - f'pick movement : {sum(w[k, m].x for k in range(K))}') + f'pick movement : {sum(z[k, m].x for k in range(K))}') pcb_part_indices = defaultdict(list) for idx, data in pcb_data.iterrows(): @@ -162,22 +171,25 @@ def line_optimizer_model(component_data, pcb_data, machine_num, hinter=True): average_pos = round( (sum(head_place_pos) / len(head_place_pos) + stopper_pos[0] - slotf1_pos[0] + 1) / slot_interval) - print(f'average_pos: {average_pos}') + for k in range(len(feeder_slot_result)): for h in range(H): if feeder_slot_result[k][h] == -1: continue feeder_slot_result[k][h] = feeder_slot_result[k][h] * 2 + average_pos - placement_result, head_sequence = greedy_placement_route_generation(partial_component_data, partial_pcb_data, - component_result, cycle_result, - feeder_slot_result, hinter=False) - print('----- Placement machine ' + str(m + 1) + ' ----- ') + placement_result, head_sequence = place_allocate_sequence_route_generation(partial_component_data, + partial_pcb_data, + component_result, cycle_result, + feeder_slot_result, hinter=False) + opt_res = OptResult(component_result, cycle_result, feeder_slot_result, placement_result, head_sequence) - info = placement_info_evaluation(partial_component_data, partial_pcb_data, opt_res, hinter=False) - optimization_assign_result(partial_component_data, partial_pcb_data, opt_res, nozzle_hinter=True, - component_hinter=True, feeder_hinter=True) - info.print() - assembly_info.append(info) - print('------------------------------ ') + info = placement_info_evaluation(partial_component_data, partial_pcb_data, opt_res, hinter=hinter) + if hinter: + print('----- Placement machine ' + str(m + 1) + ' ----- ') + optimization_assign_result(partial_component_data, partial_pcb_data, opt_res, nozzle_hinter=True, + component_hinter=True, feeder_hinter=True) + info.print() + assembly_info.append(info) + print('------------------------------ ') return assembly_info diff --git a/lineopt_reconfiguration.py b/lineopt_reconfiguration.py index c31d4c4..b72bf07 100644 --- a/lineopt_reconfiguration.py +++ b/lineopt_reconfiguration.py @@ -291,6 +291,328 @@ def evolutionary_component_assignment(pcb_data, component_data, machine_number, return min(pop_val), population[np.argmin(pop_val)] +class SpiderMonkeyOpt: + def __init__(self, pop_size, pcb_data, component_data, machine_number, estimator): + self.PcbData = pcb_data + self.ComponentData = component_data + self.Estimator = estimator + + self.PopSize = pop_size + self.LocalLimit = pop_size + self.GlobalLimit = pop_size + + self.MachineNum = machine_number + self.GroupSize = 0 + # self.Dim = sum(data.fdn for _, data in component_data.iterrows()) + machine_number + + self.CpPoints = defaultdict(int) + self.CpIndex = defaultdict(int) + self.CpNozzle = defaultdict(str) + + self.Dim = 0 + for cp_idx, data in component_data.iterrows(): + # cp_feeders[cp_idx] = data.fdn + + division_data = copy.deepcopy(data) + feeder_limit, total_points = division_data.fdn, division_data.points + surplus_points = total_points % feeder_limit + for _ in range(feeder_limit): + division_data.fdn, division_data.points = 1, math.floor(total_points / feeder_limit) + if surplus_points: + division_data.points += 1 + surplus_points -= 1 + + self.CpPoints[self.Dim], self.CpNozzle[self.Dim] = division_data.points, division_data.nz + self.CpIndex[self.Dim] = cp_idx + self.Dim += 1 + + self.Dim += machine_number + + component_list = list(range(len(self.CpPoints))) + self.GenPosition = [] + self.GenPopVal = [] + for _ in range(self.PopSize): + random.shuffle(component_list) + self.GenPosition.append(component_list.copy()) + idx, prev = 0, 0 + div = random.sample(range(len(component_list)), self.MachineNum - 1) + div.append(len(component_list)) + div.sort(reverse=False) + + for _ in range(1, self.MachineNum + 1): + self.GenPosition[-1].append(div[idx] - prev) + prev = div[idx] + idx += 1 + + self.GenPopVal.append(self.CalIndividualVal(self.GenPosition[-1])) + + self.GroupPoint = [[0 for _ in range(2)] for _ in range(self.PopSize)] + self.GroupPart = 1 + + self.GenProb = None + self.GlobalMin = self.GenPopVal[0] + self.GlobalLeaderPosition = self.GenPosition[0].copy() + self.GlobalLimitCount = 0 + + self.LocalMin = np.ones(self.PopSize) * 1e10 + self.LocalLeaderPosition = [[0 for _ in range(self.Dim)] for _ in range(self.PopSize)] + self.LocalLimitCount = [0 for _ in range(self.PopSize)] + + for k in range(self.GroupSize): + self.LocalMin[k] = self.GenPopVal[int(self.GroupPoint[k, 0])] + self.LocalLeaderPosition[k,:] = self.GenPosition[int(self.GroupPoint[k, 0]),:] + + self.CrossoverRatio = 0.1 + + + def GlobalLearning(self): + GlobalTrial = self.GlobalMin + for i in range(self.PopSize): + if self.GenPopVal[i] < self.GlobalMin: + self.GlobalMin = self.GenPopVal[i] + self.GlobalLeaderPosition = self.GenPosition[i].copy() + + if math.fabs(GlobalTrial - self.GlobalMin) < 1e-5: + self.GlobalLimitCount = self.GlobalLimitCount + 1 + else: + self.GlobalLimitCount = 0 + + def LocalLearning(self): + OldMin = np.zeros(self.PopSize) + for k in range(self.GroupSize): + OldMin[k] = self.LocalMin[k] + + for k in range(self.GroupSize): + i = int(self.GroupPoint[k][0]) + while i <= int(self.GroupPoint[k][1]): + if self.GenPopVal[i] < self.LocalMin[k]: + self.LocalMin[k] = self.GenPopVal[i] + self.LocalLeaderPosition[k] = self.GenPosition[i].copy() + i = i + 1 + + for k in range(self.GroupSize): + if math.fabs(OldMin[k] - self.LocalMin[k]) < 1e-5: + self.LocalLimitCount[k] = self.LocalLimitCount[k] + 1 + else: + self.LocalLimitCount[k] = 0 + + def CalculateProbabilities(self): + self.GenProb = [0 for _ in range(self.PopSize)] + MaxVal = self.GenPopVal[0] + i = 1 + while i < self.PopSize: + if self.GenPopVal[i] > MaxVal: + MaxVal = self.GenPopVal[i] + i += 1 + for i in range(self.PopSize): + self.GenProb[i] = (0.9 * (self.GenPopVal[i] / MaxVal)) + 0.1 + + def LocalLeaderPhase(self, k): + lo = int(self.GroupPoint[k][0]) + hi = int(self.GroupPoint[k][1]) + i = lo + while i <= hi: + + NewGene1, NewGene2 = self.CrossoverOperation(self.GenPosition[i], self.LocalLeaderPosition[k]) + NewGeneVal1, NewGeneVal2 = self.CalIndividualVal(NewGene1), self.CalIndividualVal(NewGene2) + + if NewGeneVal1 < self.GenPopVal[i]: + self.GenPosition[i] = NewGene1 + self.GenPopVal[i] = NewGeneVal1 + + if NewGeneVal2 < self.GenPopVal[i]: + self.GenPosition[i] = NewGene2 + self.GenPopVal[i] = NewGeneVal2 + i += 1 + + def GlobalLeaderPhase(self, k): + lo = int(self.GroupPoint[k][0]) + hi = int(self.GroupPoint[k][1]) + i = lo + l = lo + while l < hi: + if random.random() < self.GenProb[i]: + l += 1 + NewGene1, NewGene2 = self.CrossoverOperation(self.GenPosition[i], self.GlobalLeaderPosition) + NewGeneVal1, NewGeneVal2 = self.CalIndividualVal(NewGene1), self.CalIndividualVal(NewGene2) + + if NewGeneVal1 < self.GenPopVal[i]: + self.GenPosition[i] = NewGene1 + self.GenPopVal[i] = NewGeneVal1 + + if NewGeneVal2 < self.GenPopVal[i]: + self.GenPosition[i] = NewGene2 + self.GenPopVal[i] = NewGeneVal2 + + i += 1 + if i == hi: + i = lo + + def LocalLeaderDecision(self): + for k in range(self.GroupSize): + if self.LocalLimitCount[k] > self.LocalLimit: + i = self.GroupPoint[k][0] + while i <= int(self.GroupPoint[k][1]): + if random.random() >= self.CrossoverRatio: + NewGenPosition = list(range(self.Dim - self.MachineNum)) + random.shuffle(NewGenPosition) + + idx, prev = 0, 0 + div = random.sample(range(len(NewGenPosition)), self.MachineNum - 1) + div.append(len(NewGenPosition)) + div.sort(reverse=False) + + for _ in range(1, self.MachineNum + 1): + NewGenPosition.append(div[idx] - prev) + prev = div[idx] + idx += 1 + + NewGenVal = self.CalIndividualVal(NewGenPosition) + + if NewGenVal < self.GenPopVal[i]: + self.GenPosition[i] = NewGenPosition.copy() + self.GenPopVal[i] = NewGenVal + else: + NewGene1, NewGene2 = self.CrossoverOperation(self.GenPosition[i], self.GlobalLeaderPosition) + NewGeneVal1, NewGeneVal2 = self.CalIndividualVal(NewGene1), self.CalIndividualVal(NewGene2) + + if NewGeneVal1 < self.GenPopVal[i]: + self.GenPosition[i] = NewGene1.copy() + self.GenPopVal[i] = NewGeneVal1 + + if NewGeneVal2 < self.GenPopVal[i]: + self.GenPosition[i] = NewGene2.copy() + self.GenPopVal[i] = NewGeneVal2 + + i += 1 + + self.LocalLimitCount[k] = 0 + + def GlobalLeaderDecision(self): + if self.GlobalLimitCount> self.GlobalLimit: + self.GroupPart += 1 + self.GlobalLimitCount = 0 + self.CreateGroup() + self.LocalLearning() + + def CreateGroup(self): + g = 0 + lo = 0 + + while lo < self.PopSize: + hi = lo + int(self.PopSize / self.GroupPart) + self.GroupPoint[g][0] = lo + self.GroupPoint[g][1] = hi + if self.PopSize - hi < int(self.PopSize / self.GroupPart): + self.GroupPoint[g][1] = (self.PopSize - 1) + g = g + 1 + lo = hi + 1 + self.GroupSize = g + + def CrossoverOperation(self, gene1, gene2): + len_ = len(gene1) + sub1, sub2 = partially_mapped_crossover(gene1[0: len_ - self.MachineNum], gene2[0: len_ - self.MachineNum]) + + pos1, pos2 = random.randint(0, self.MachineNum - 1), random.randint(0, self.MachineNum - 1) + machine_assign1, machine_assign2 = gene1[len_ - self.MachineNum:], gene2[len_ - self.MachineNum:] + machine_assign1[pos1], machine_assign2[pos2] = machine_assign2[pos2], machine_assign1[pos1] + while sum(machine_assign1) != len_ - self.MachineNum: + machine_idx = random.randint(0, self.MachineNum - 1) + if machine_assign1[machine_idx] == 0: + continue + if sum(machine_assign1) > len_ - self.MachineNum: + machine_assign1[machine_idx] -= 1 + else: + machine_assign1[machine_idx] += 1 + + while sum(machine_assign2) != len_ - self.MachineNum: + machine_idx = random.randint(0, self.MachineNum - 1) + if machine_assign2[machine_idx] == 0: + continue + if sum(machine_assign2) > len_ - self.MachineNum: + machine_assign2[machine_idx] -= 1 + else: + machine_assign2[machine_idx] += 1 + + sub1.extend(machine_assign1) + sub2.extend(machine_assign2) + return sub1, sub2 + + + def CalIndividualVal(self, gene): + ComponentNum = len(self.ComponentData) + assignment_result = [[0 for _ in range(ComponentNum)] for _ in range(self.MachineNum)] + + idx, machine_index = 0, 0 + for num in range(self.Dim - self.MachineNum, self.Dim): + for _ in range(gene[num]): + assignment_result[machine_index][self.CpIndex[gene[idx]]] += self.CpPoints[gene[idx]] + idx += 1 + machine_index += 1 + + val = 0 + cp_items = converter(self.PcbData, self.ComponentData, assignment_result) + for machine_index in range(self.MachineNum): + cp_points, cp_nozzle, board_width, board_height = cp_items[machine_index] + val = max(val, self.Estimator.predict(cp_points, cp_nozzle, board_width, board_height)) + return val + +def spider_monkey_component_assignment(pcb_data, component_data, machine_number, estimator): + population_size, iteration_number = 20, 50 + + smo = SpiderMonkeyOpt(population_size, pcb_data, component_data, machine_number, estimator) + + # ========================== Calling: GlobalLearning() ======================== # + smo.GlobalLearning() + + # ========================== Calling: create_group() ========================== # + smo.CreateGroup() + + # ========================= Calling: LocalLearning() ========================== # + smo.LocalLearning() + + # ================================= Looping ================================== # + with tqdm(total=iteration_number) as pbar: + pbar.set_description('spider monkey algorithm process for PCB assembly line balance') + for _ in range(iteration_number): + for k in range(smo.GroupSize): + # ==================== Calling: LocalLeaderPhase() =================== # + smo.LocalLeaderPhase(k) + + # =================== Calling: CalculateProbabilities() ================== # + smo.CalculateProbabilities() + + for k in range(smo.GroupSize): + # ==================== Calling: GlobalLeaderPhase() ================== # + smo.GlobalLeaderPhase(k) + + # ======================= Calling: GlobalLearning() ====================== # + smo.GlobalLearning() + + # ======================= Calling: LocalLearning() ======================= # + smo.LocalLearning() + + # ================== Calling: LocalLeaderDecision() ====================== # + smo.LocalLeaderDecision() + + # ===================== Calling: GlobalLeaderDecision() ================== # + smo.GlobalLeaderDecision() + + pbar.update(1) + + assignment_result = [[0 for _ in range(len(component_data))] for _ in range(machine_number)] + + idx, machine_index = 0, 0 + for num in range(smo.Dim - machine_number, smo.Dim): + for _ in range(smo.GlobalLeaderPosition[num]): + assignment_result[machine_index][smo.CpIndex[smo.GlobalLeaderPosition[idx]]] += \ + smo.CpPoints[smo.GlobalLeaderPosition[idx]] + idx += 1 + machine_index += 1 + + return smo.GlobalMin, assignment_result + + @timer_wrapper def line_optimizer_reconfiguration(component_data, pcb_data, machine_number): # === assignment of heads to modules is omitted === @@ -303,9 +625,8 @@ def line_optimizer_reconfiguration(component_data, pcb_data, machine_number): print('random component allocation algorithm process for PCB assembly line balance') val, assignment = random_component_assignment(pcb_data, component_data, machine_number, estimator) elif i == 1: - # brute force - # which is proved to be useless, since it only ran in reasonable time for the smaller test instances - continue + # spider monkey + val, assignment = spider_monkey_component_assignment(pcb_data, component_data, machine_number, estimator) elif i == 2: # local search print('local search component allocation algorithm process for PCB assembly line balance') @@ -314,7 +635,8 @@ def line_optimizer_reconfiguration(component_data, pcb_data, machine_number): # evolutionary val, assignment = evolutionary_component_assignment(pcb_data, component_data, machine_number, estimator) else: - # greedy: unclear description + # brute force + # which is proved to be useless, since it only ran in reasonable time for the smaller test instances continue if optimal_val is None or val < optimal_val: @@ -324,3 +646,5 @@ def line_optimizer_reconfiguration(component_data, pcb_data, machine_number): raise Exception('no feasible solution! ') return optimal_assignment + + diff --git a/optimizer.py b/optimizer.py index 225af6e..dee9b40 100644 --- a/optimizer.py +++ b/optimizer.py @@ -1,7 +1,3 @@ -import time - -import pandas as pd - from dataloader import * from lineopt_genetic import line_optimizer_genetic from lineopt_heuristic import line_optimizer_heuristic @@ -12,9 +8,9 @@ from lineopt_model import line_optimizer_model from base_optimizer.optimizer_interface import * -def optimizer(pcb_data, component_data, feeder_data, params): +def optimizer(pcb_data, component_data, feeder_data, params, hinter=True): if params.machine_number == 1: - assembly_info = [base_optimizer(1, pcb_data, component_data, feeder_data, params, hinter=True)] + assembly_info = [base_optimizer(1, pcb_data, component_data, feeder_data, params, hinter=hinter)] return assembly_info if params.line_optimizer == 'hyper-heuristic' or params.line_optimizer == 'heuristic' or params.line_optimizer \ @@ -33,7 +29,7 @@ def optimizer(pcb_data, component_data, feeder_data, params): for machine_index in range(params.machine_number): assembly_info.append(base_optimizer(machine_index + 1, partial_pcb_data[machine_index], partial_component_data[machine_index], - pd.DataFrame(columns=['slot', 'part']), params, hinter=True)) + pd.DataFrame(columns=['slot', 'part']), params, hinter=hinter)) elif params.line_optimizer == 'mip-model': assembly_info = line_optimizer_model(component_data, pcb_data, params.machine_number) else: @@ -49,16 +45,14 @@ def main(): parser.add_argument('--mode', default=1, type=int, help='mode: 0 -directly load pcb data without optimization ' 'for data analysis, 1 -optimize pcb data, 2 -batch test') parser.add_argument('--filename', default='PCB.txt', type=str, help='load pcb data') - # parser.add_argument('--filename', default='chapter3-2/PCB2-8 Arg1.txt', type=str, help='load pcb data') + parser.add_argument('--comp_register', default=1, type=int, help='register the component according the pcb data') - parser.add_argument('--machine_number', default=1, type=int, help='the number of machine in the assembly line') - # parser.add_argument('--machine_optimizer', default='mip-model', type=str, help='optimizer for single machine') + parser.add_argument('--machine_number', default=3, type=int, help='the number of machine in the assembly line') parser.add_argument('--machine_optimizer', default='feeder-priority', type=str, help='optimizer for single machine') parser.add_argument('--line_optimizer', default='hyper-heuristic', type=str, help='optimizer for PCB assembly line') - # parser.add_argument('--line_optimizer', default='model', type=str, help='optimizer for PCB assembly line') parser.add_argument('--feeder_limit', default=1, type=int, help='the upper feeder limit for each type of component') parser.add_argument('--save', default=0, type=int, help='save the optimization result') - parser.add_argument('--save_suffix', default='(10)', type=str, help='load pcb data') + parser.add_argument('--save_suffix', default='(1)', type=str, help='load pcb data') params = parser.parse_args() # 结果输出显示所有行和列 @@ -69,6 +63,7 @@ def main(): assembly_info = [] for machine_index in range(len(partial_pcb_data)): opt_res = convert_pcbdata_to_result(partial_pcb_data[machine_index], partial_component_data[machine_index]) + info = placement_info_evaluation(partial_component_data[machine_index], partial_pcb_data[machine_index], opt_res) assembly_info.append(info) @@ -76,6 +71,7 @@ def main(): nozzle_hinter=True, component_hinter=True, feeder_hinter=True) info.print() + if params.save: output_optimize_result(f'result/{params.filename[:-4]}-T-Solution-M0{machine_index + 1}', partial_component_data[machine_index], partial_pcb_data[machine_index], opt_res) @@ -89,7 +85,6 @@ def main(): f'standard deviation: {np.std([info.total_time for info in assembly_info]): .3f}') elif params.mode == 1: - # sys.stdout = open(f'record/{params.filename[:-4]}-{params.line_optimizer}.txt', 'w') # 加载PCB数据 partial_pcb_data, partial_component_data, feeder_data = load_data(params.filename, load_feeder=True) @@ -108,28 +103,33 @@ def main(): print(f'finial assembly time: {max(info.total_time for info in assembly_info): .3f} s, ' f'standard deviation: {np.std([info.total_time for info in assembly_info]): .3f}') + elif params.mode == 2: + + # sys.stdout = open(f'record/dissertation-experiment.txt', 'w') machine_optimizer = ['two-phase', 'hybrid-genetic', 'cell-division', 'feeder-priority', 'aggregation'] running_round = 10 opt_columns = ['Cycle', 'Pick', 'Nozzle-Change', 'Running-Time'] opt_result, opt_runtime = defaultdict(pd.DataFrame), defaultdict(pd.DataFrame) for opt in machine_optimizer: - opt_result[opt] =pd.DataFrame(columns=opt_columns) + opt_result[opt] = pd.DataFrame(columns=opt_columns) opt_result[opt].index.name = 'file' for _, file in enumerate(os.listdir('data/')): if file[-3:] != 'txt': continue + partial_pcb_data, partial_component_data, feeder_data = load_data(file) pcb_data, component_data = merge_data(partial_pcb_data, partial_component_data) + for opt in machine_optimizer: for round_idx in range(running_round): - print(f'--- file ： {file}, round ： {round_idx}, optimizer ： {opt} --- ') + print(f'--- file ： {file}, round ： {round_idx + 1}, optimizer ： {opt} --- ') - params = parser.parse_args(['--machine_optimizer', opt, '--machine_number', str(1)]) + params = parser.parse_args(['--machine_optimizer', opt, '--machine_number', str(1), '--filename', file]) start_time = time.time() - assembly_info = optimizer(pcb_data, component_data, feeder_data, params) + assembly_info = optimizer(pcb_data, component_data, feeder_data, params, hinter=False) opt_result[opt].loc[file + str(round_idx + 1), 'Cycle'] = assembly_info[0].cycle_counter opt_result[opt].loc[file + str(round_idx + 1), 'Pick'] = assembly_info[0].pickup_counter @@ -140,8 +140,7 @@ def main(): for opt, result in opt_result.items(): result.to_excel(writer, sheet_name=opt, float_format='%.3f', na_rep='') else: - # line_optimizer = ['T-Solution', 'hyper-heuristic', 'genetic', 'reconfiguration'] - line_optimizer = ['genetic'] + line_optimizer = ['T-Solution', 'hyper-heuristic', 'genetic', 'reconfiguration'] file_dirs = ['L01', 'L02', 'L03'] running_round = 10 @@ -169,6 +168,7 @@ def main(): warning_info = f'file: {file_dir}/{file}: an unexpected error occurs for data loader' warnings.warn(warning_info, SyntaxWarning) continue + machine_number = len(partial_pcb_data) if not os.path.exists(f'record/{file_dir}'): os.makedirs(f'record/{file_dir}') @@ -207,19 +207,20 @@ def main(): ['--filename', file_dir + '/' + file, '--machine_number', str(machine_number), '--line_optimizer', line_opt, '--save_suffix', f'({round_idx + 1})']) start_time = time.time() - assembly_info = optimizer(merge_pcb_data, merge_component_data, params) + assembly_info = optimizer(merge_pcb_data, merge_component_data, None, params) line_opt_result[file_dir].loc[file, line_opt + str(round_idx + 1)] = max( info.total_time for info in assembly_info) line_opt_runtime[file_dir].loc[file, line_opt + str(round_idx + 1)] = time.time() - start_time - for machine_idx, info in enumerate(assembly_info): - print( - f'assembly time for machine {machine_idx + 1: d}: {info.total_time: .3f} s, ' - f'total placement: {info.total_points}, ' - f'total component types {info.total_components: d}') - print(f'finial assembly time: {max(info.total_time for info in assembly_info): .3f} s, ' - f'standard deviation: {np.std([info.total_time for info in assembly_info]): .3f}') + for machine_idx, info in enumerate(assembly_info): + print( + f'assembly time for machine {machine_idx + 1: d}: {info.total_time: .3f} s, ' + f'total placement: {info.total_points}, ' + f'total component types {info.total_components: d}') + + print(f'finial assembly time: {max(info.total_time for info in assembly_info): .3f} s, ' + f'standard deviation: {np.std([info.total_time for info in assembly_info]): .3f}') with pd.ExcelWriter('result/line_optimizer.xlsx', engine='openpyxl') as writer: for file_dir, result in line_opt_result.items():