调整工程架构，增补了几种算法，初步添加神经网路训练拟合代码

.gitignore

@@ -6,3 +6,6 @@ __pycache__
 Lib/
 Scripts/
 *.cfg
+*.pkl
+*.pth
+*.m

base_optimizer/optimizer_aggregation.py

@@ -1,8 +1,5 @@
 from base_optimizer.optimizer_common import *
 
-from gurobipy import *
-from collections import defaultdict
-
 
 def list_range(start, end=None):
     return list(range(start)) if end is None else list(range(start, end))

base_optimizer/optimizer_celldivision.py

@@ -1,5 +1,30 @@
 from base_optimizer.optimizer_common import *
-from result_analysis import *
+
+
+def component_assign_evaluate(component_data, component_result, cycle_result, feeder_slot_result) -> float:
+    nozzle_change_counter = 0
+    for head in range(max_head_index):
+        nozzle = ''
+        for cycle in range(len(component_result)):
+            component_index = component_result[cycle][head]
+            if component_index == -1:
+                continue
+
+            if cycle != 0 and nozzle != component_data.loc[component_index, 'nz']:
+                nozzle_change_counter += 1
+            nozzle = component_data.loc[component_index, 'nz']
+
+    gang_pick_counter = 0
+    for cycle, feeder_slot in enumerate(feeder_slot_result):
+        pick_slot = defaultdict(int)
+        for head, slot in enumerate(feeder_slot):
+            if slot == -1:
+                continue
+            pick_slot[slot - head * interval_ratio] += 1
+        for _ in pick_slot.values():
+            gang_pick_counter += cycle_result[cycle]
+
+    return sum(cycle_result) + e_nz_change * nozzle_change_counter + e_gang_pick * gang_pick_counter
 
 
 def convert_cell_2_result(pcb_data, component_data, component_cell, population):

base_optimizer/optimizer_common.py

@@ -1,25 +1,22 @@
-import copy
-import time
-import math
-import random
-import argparse
-import os
-import warnings
-import copy
-
-import numpy as np
-import pandas as pd
-import matplotlib.pyplot as plt
-
 from functools import wraps
 from collections import defaultdict
 from tqdm import tqdm
+from gurobipy import *
+from sklearn.linear_model import LinearRegression
 
-# 整线参数
-max_machine_index = 3
-
-# 时间参数
-T_pp, T_tr, T_nc, T_pl = 2, 5, 25, 0
+import os
+import time
+import math
+import random
+import copy
+import torch
+import argparse
+import joblib
+import pickle
+import warnings
+import numpy as np
+import pandas as pd
+import matplotlib.pyplot as plt
 
 # 机器参数
 max_head_index, max_slot_index = 6, 120
@@ -56,6 +53,73 @@ t_nozzle_put, t_nozzle_pick = 0.9, 0.75  # 装卸吸嘴用时
 t_nozzle_change = t_nozzle_put + t_nozzle_pick
 t_fix_camera_check = 0.12  # 固定相机检测时间
 
+# 时间参数（整线相关）
+T_pp, T_tr, T_nc, T_pl = 2, 5, 25, 0
+
+
+class OptInfo:
+    def __init__(self):
+        self.placement_time = 0
+
+        self.cycle_counter = 0
+        self.nozzle_change_counter = 0
+        self.pickup_counter = 0
+
+        self.pickup_movement = 0
+        self.placement_movement = 0
+
+
+def optimization_assign_result(component_data, pcb_data, component_result, cycle_result, feeder_slot_result,
+                               nozzle_hinter=False, component_hinter=False, feeder_hinter=False):
+    if nozzle_hinter:
+        columns = ['H{}'.format(i + 1) for i in range(max_head_index)] + ['cycle']
+
+        nozzle_assign = pd.DataFrame(columns=columns)
+        for cycle, components in enumerate(component_result):
+            nozzle_assign.loc[cycle, 'cycle'] = cycle_result[cycle]
+            for head in range(max_head_index):
+                index = component_result[cycle][head]
+                if index == -1:
+                    nozzle_assign.loc[cycle, 'H{}'.format(head + 1)] = ''
+                else:
+                    nozzle_assign.loc[cycle, 'H{}'.format(head + 1)] = component_data.loc[index].nz
+
+        print(nozzle_assign)
+        print('')
+
+    if component_hinter:
+        columns = ['H{}'.format(i + 1) for i in range(max_head_index)] + ['cycle']
+
+        component_assign = pd.DataFrame(columns=columns)
+        for cycle, components in enumerate(component_result):
+            component_assign.loc[cycle, 'cycle'] = cycle_result[cycle]
+            for head in range(max_head_index):
+                index = component_result[cycle][head]
+                if index == -1:
+                    component_assign.loc[cycle, 'H{}'.format(head + 1)] = ''
+                else:
+                    component_assign.loc[cycle, 'H{}'.format(head + 1)] = component_data.loc[index].part
+
+        print(component_assign)
+        print('')
+
+    if feeder_hinter:
+        columns = ['H{}'.format(i + 1) for i in range(max_head_index)] + ['cycle']
+
+        feedr_assign = pd.DataFrame(columns=columns)
+        for cycle, components in enumerate(feeder_slot_result):
+            feedr_assign.loc[cycle, 'cycle'] = cycle_result[cycle]
+            for head in range(max_head_index):
+                slot = feeder_slot_result[cycle][head]
+                if slot == -1:
+                    feedr_assign.loc[cycle, 'H{}'.format(head + 1)] = 'A'
+                else:
+                    feedr_assign.loc[cycle, 'H{}'.format(head + 1)] = 'F{}'.format(
+                        slot) if slot <= max_slot_index // 2 else 'R{}'.format(slot - max_head_index)
+
+        print(feedr_assign)
+        print('')
+
 
 def axis_moving_time(distance, axis=0):
     distance = abs(distance) * 1e-3
@@ -315,8 +379,8 @@ def dynamic_programming_cycle_path(pcb_data, cycle_placement, assigned_feeder):
             print(assigned_feeder)
             print(cycle_placement)
 
-        pos.append([pcb_data.loc[placement]['x'] - head * head_interval + stopper_pos[0],
-                    pcb_data.loc[placement]['y'] + stopper_pos[1]])
+        pos.append([pcb_data.iloc[placement]['x'] - head * head_interval + stopper_pos[0],
+                    pcb_data.iloc[placement]['y'] + stopper_pos[1]])
 
         feeder_set.add(feeder - head * interval_ratio)
 
@@ -369,8 +433,8 @@ def dynamic_programming_cycle_path(pcb_data, cycle_placement, assigned_feeder):
         head_sequence.append(head_set[parent - 1])
 
     start_head, end_head = head_sequence[0], head_sequence[-1]
-    if pcb_data.loc[cycle_placement[start_head]]['x'] - start_head * head_interval > \
-            pcb_data.loc[cycle_placement[end_head]]['x'] - end_head * head_interval:
+    if pcb_data.iloc[cycle_placement[start_head]]['x'] - start_head * head_interval > \
+            pcb_data.iloc[cycle_placement[end_head]]['x'] - end_head * head_interval:
         head_sequence = list(reversed(head_sequence))
     return head_sequence
 
@@ -382,11 +446,11 @@ def greedy_placement_route_generation(component_data, pcb_data, component_result
     mount_point_pos = [[] for _ in range(len(component_data))]
 
     for i in range(len(pcb_data)):
-        part = pcb_data.loc[i]['part']
+        part = pcb_data.iloc[i]['part']
         component_index = component_data[component_data['part'] == part].index.tolist()[0]
         # 记录贴装点序号索引和对应的位置坐标
         mount_point_index[component_index].append(i)
-        mount_point_pos[component_index].append([pcb_data.loc[i]['x'], pcb_data.loc[i]['y']])
+        mount_point_pos[component_index].append([pcb_data.iloc[i]['x'], pcb_data.iloc[i]['y']])
 
     search_dir = 1  # 0：自左向右搜索  1：自右向左搜索
     for cycle_set in range(len(component_result)):
@@ -883,18 +947,18 @@ def swap_mutation(parent):
     return parent
 
 
-def constraint_swap_mutation(component_points, individual):
+def constraint_swap_mutation(component_points, individual, machine_number):
     offspring = individual.copy()
 
     idx, component_index = 0, random.randint(0, len(component_points) - 1)
     for _, points in component_points:
         if component_index == 0:
             while True:
-                index1, index2 = random.sample(range(points + max_machine_index - 2), 2)
+                index1, index2 = random.sample(range(points + machine_number - 2), 2)
                 if offspring[idx + index1] != offspring[idx + index2]:
                     break
 
-            clip = offspring[idx: idx + points + max_machine_index - 1].copy()
+            clip = offspring[idx: idx + points + machine_number - 1].copy()
             avl_machine = 0
             for idx_, gene in enumerate(clip):
                 if gene == 0 and (idx_ == 0 or clip[idx_ - 1] != 0):
@@ -912,7 +976,7 @@ def constraint_swap_mutation(component_points, individual):
             break
 
         component_index -= 1
-        idx += (points + max_machine_index - 1)
+        idx += (points + machine_number - 1)
 
     return offspring
 
@@ -959,3 +1023,92 @@ def get_top_k_value(pop_val, k: int, reverse=True):
                 res.append(j)
                 break
     return res
+
+
+def get_line_config_number(machine_number, component_number):
+    div_counter = 0
+    div_set = set()
+    for div1 in range(component_number - 2):
+        for div2 in range(div1 + 1, component_number - 1):
+            machine_div = [div1 + 1, div2 - div1, component_number - div2 - 1]
+            machine_div.sort()
+            div_str = "".join(str(s) + '|' for s in machine_div)
+            if div_str in div_set:
+                continue
+
+            div_set.add(div_str)
+            assign_component_counter = defaultdict(list)
+            for div in machine_div:
+                assign_component_counter[div] += 1
+
+            case_div_counter, case_comp_number = 1, component_number
+            for idx in range(machine_number - 1):
+                div = 1
+                while machine_div[idx]:
+                    div *= (case_comp_number / machine_div[idx])
+                    case_comp_number -= 1
+                    machine_div[idx] -= 1
+
+                case_div_counter *= div
+
+            for key, val in assign_component_counter.values():
+                div = 1
+                while val:
+                    div *= val
+                    val -= 1
+                case_div_counter /= div
+            div_counter += case_div_counter
+        return div_counter
+
+
+def partial_data_convert(pcb_data, component_data, machine_assign, machine_number):
+    assignment_result = copy.deepcopy(machine_assign)
+    partial_pcb_data, partial_component_data = defaultdict(pd.DataFrame), defaultdict(pd.DataFrame)
+    for machine_index in range(machine_number):
+        partial_pcb_data[machine_index] = pd.DataFrame(columns=pcb_data.columns)
+        partial_component_data[machine_index] = component_data.copy(deep=True)
+
+    # === averagely assign available feeder ===
+    for part_index, data in component_data.iterrows():
+        feeder_limit = data['feeder-limit']
+        feeder_points = [assignment_result[machine_index][part_index] for machine_index in range(max_machine_index)]
+
+        for machine_index in range(machine_number):
+            if feeder_points[machine_index] == 0:
+                continue
+
+            arg_feeder = max(math.floor(feeder_points[machine_index] / sum(feeder_points) * data['feeder-limit']), 1)
+
+            partial_component_data[machine_index].loc[part_index]['feeder-limit'] = arg_feeder
+            feeder_limit -= arg_feeder
+
+        for machine_index in range(machine_number):
+            if feeder_limit <= 0:
+                break
+
+            if feeder_points[machine_index] == 0:
+                continue
+            partial_component_data[machine_index].loc[part_index]['feeder-limit'] += 1
+            feeder_limit -= 1
+
+        for machine_index in range(machine_number):
+            if feeder_points[machine_index] > 0:
+                assert partial_component_data[machine_index].loc[part_index]['feeder-limit'] > 0
+
+    # === assign placements ===
+    component_machine_index = [0 for _ in range(len(component_data))]
+    for _, data in pcb_data.iterrows():
+        part_index = component_data[component_data['part'] == data['part']].index.tolist()[0]
+        while True:
+            machine_index = component_machine_index[part_index]
+            if assignment_result[machine_index][part_index] == 0:
+                component_machine_index[part_index] += 1
+                machine_index += 1
+            else:
+                break
+        assignment_result[machine_index][part_index] -= 1
+        partial_pcb_data[machine_index] = pd.concat([partial_pcb_data[machine_index], pd.DataFrame(data).T])
+
+    return partial_pcb_data, partial_component_data
+
+

base_optimizer/optimizer_feederpriority.py

@@ -13,20 +13,20 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
     feeder_base = [-2] * max_slot_index   # 已安装在供料器基座上的元件（-2: 未分配，-1: 占用状态）
     feeder_base_points = [0] * max_slot_index   # 供料器基座结余贴装点数量
 
-    for data in pcb_data.iterrows():
-        pos, part = data[1]['x'] + stopper_pos[0], data[1]['part']
+    for _, data in pcb_data.iterrows():
+        pos, part = data.x + stopper_pos[0], data.part
 
-        part_index = component_data[component_data['part'] == part].index.tolist()[0]
+        part_index = component_data[component_data.part == part].index.tolist()[0]
         if part not in component_data:
             feeder_limit[part_index] = component_data.loc[part_index]['feeder-limit']
             feeder_arrange[part_index] = 0
 
         feeder_points[part_index] += 1
         mount_center_pos[part_index] += ((pos - mount_center_pos[part_index]) / feeder_points[part_index])
-        part_nozzle[part_index] = component_data.loc[part_index]['nz']
+        part_nozzle[part_index] = component_data.loc[part_index].nz
 
     for part_index, points in feeder_points.items():
-        feeder_division_points[part_index] = max(points // feeder_limit[part_index], 1)
+        feeder_division_points[part_index] = points // feeder_limit[part_index]
 
     nozzle_component, nozzle_component_points = defaultdict(list), defaultdict(list)
     for part, nozzle in part_nozzle.items():
@@ -36,13 +36,13 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
 
     if feeder_data is not None:
         for _, feeder in feeder_data.iterrows():
-            slot, part = feeder['slot'], feeder['part']
-            part_index = component_data[component_data['part'] == part].index.tolist()[0]
+            slot, part = feeder.slot, feeder.part
+            part_index = component_data[component_data.part == part].index.tolist()[0]
 
             # 供料器基座分配位置和对应贴装点数
             feeder_base[slot], feeder_base_points[slot] = part_index, feeder_division_points[part_index]
 
-            feeder_type = component_data.loc[part_index]['fdr']
+            feeder_type = component_data.loc[part_index].fdr
             extra_width = feeder_width[feeder_type][0] + feeder_width[feeder_type][1] - slot_interval
             while extra_width > 0:
                 slot += 1
@@ -75,15 +75,60 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
     nozzle_pattern_list.sort(key=lambda x: x[1], reverse=True)
 
     # 后确定吸嘴分配模式
-    head_index = [3, 2, 4, 1, 5, 0]
-    nozzle_pattern = [0] * max_head_index
-    for nozzle, _ in nozzle_pattern_list:
-        counter = nozzle_assigned_counter[nozzle]
+    upper_head, extra_head = defaultdict(int), defaultdict(int)
+    head_index = []
+    for nozzle, head in nozzle_assigned_counter.items():
+        # 每个吸嘴能达成同时拾取数目的上限
+        upper_head[nozzle] = min(len(nozzle_component[nozzle]), head)
+        extra_head[nozzle] = head - upper_head[nozzle]
+
+    head_counter = (sum(upper_head.values()) - 1) // 2
+    while head_counter >= 0:
+        if head_counter != (sum(upper_head.values()) - 1) - head_counter:
+            head_index.append((sum(upper_head.values()) - 1) - head_counter)
+        head_index.append(head_counter)
+        head_counter -= 1
+
+    nozzle_pattern = [None for _ in range(sum(upper_head.values()))]
+    for nozzle in upper_head.keys():
+        counter = upper_head[nozzle]
         while counter:
             nozzle_pattern[head_index[0]] = nozzle
             counter -= 1
             head_index.pop(0)
 
+    head = 0
+    while head + sum(extra_head.values()) <= len(nozzle_pattern):
+        extra_head_cpy = copy.deepcopy(extra_head)
+        increment = 0
+        while increment < sum(extra_head.values()):
+            extra_head_cpy[nozzle_pattern[head + increment]] -= 1
+            increment += 1
+
+        check_extra_head = True
+        for head_ in extra_head_cpy.values():
+            if head_ != 0:
+                check_extra_head = False    # 任一项不为0， 说明不构成
+                break
+
+        if check_extra_head:
+            increment = 0
+            while increment < sum(extra_head.values()):
+                nozzle_pattern.append(nozzle_pattern[head + increment])
+                increment += 1
+
+            for nozzle in extra_head.keys():
+                extra_head[nozzle] = 0
+
+            break
+        head += 1
+
+    for nozzle, head_ in extra_head.items():
+        while head_:
+            nozzle_pattern.append(nozzle)
+            head_ -= 1
+
+    assert len(nozzle_pattern) == max_head_index
     while True:
         best_assign, best_assign_points = [], []
         best_assign_slot, best_assign_value = -1, -np.Inf
@@ -162,7 +207,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
 
                     part = tmp_nozzle_component[nozzle_assign][index_]
 
-                    feeder_type = component_data.loc[part]['fdr']
+                    feeder_type = component_data.loc[part].fdr
                     extra_width, extra_slot = feeder_width[feeder_type][0] + feeder_width[feeder_type][1] - slot_interval, 1
                     slot_overlap = False
                     while extra_width > 0:
@@ -210,7 +255,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
                         continue
                     for idx, part in enumerate(assign_part_stack):
 
-                        feeder_type = component_data.loc[part]['fdr']
+                        feeder_type = component_data.loc[part].fdr
                         extra_width, extra_slot = feeder_width[feeder_type][0] + feeder_width[feeder_type][
                             1] - slot_interval, 1
 
@@ -223,7 +268,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
                             extra_width -= slot_interval
                             extra_slot += 1
 
-                        if component_data.loc[part]['nz'] == nozzle_pattern[head] and not slot_overlap:
+                        if component_data.loc[part].nz == nozzle_pattern[head] and not slot_overlap:
                             feeder_assign[head], feeder_assign_points[head] = assign_part_stack[idx], \
                                                                               assign_part_stack_points[idx]
 
@@ -237,7 +282,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
                     continue
                 part, points = assign_part_stack[0], assign_part_stack_points[0]
 
-                feeder_type = component_data.loc[part]['fdr']
+                feeder_type = component_data.loc[part].fdr
                 extra_width, extra_slot = feeder_width[feeder_type][0] + feeder_width[feeder_type][1] - slot_interval, 1
 
                 slot_overlap = False
@@ -259,7 +304,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
                         extra_width -= head_interval
                 else:
                     # 返还由于机械限位无法分配的，压入元件堆栈中的元素
-                    nozzle = component_data.loc[part]['nz']
+                    nozzle = component_data.loc[part].nz
                     tmp_nozzle_component[nozzle].insert(0, part)
                     tmp_nozzle_component_points[nozzle].insert(0, points)
 
@@ -269,7 +314,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
             # 仍然存在由于机械限位，无法进行分配的在堆栈中的元件
             while assign_part_stack:
                 part, points = assign_part_stack[0], assign_part_stack_points[0]
-                nozzle = component_data.loc[part]['nz']
+                nozzle = component_data.loc[part].nz
 
                 tmp_nozzle_component[nozzle].insert(0, part)
                 tmp_nozzle_component_points[nozzle].insert(0, points)
@@ -283,7 +328,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
                     continue
                 average_slot.append(
                     (mount_center_pos[feeder_] - slotf1_pos[0]) / slot_interval + 1 - head * interval_ratio)
-                if nozzle_pattern and component_data.loc[feeder_]['nz'] != nozzle_pattern[head]:
+                if nozzle_pattern and component_data.loc[feeder_].nz != nozzle_pattern[head]:
                     nozzle_change_counter += 1
 
             if len(average_slot) == 0:
@@ -349,7 +394,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
             # 更新供料器基座信息
             feeder_base[best_assign_slot + idx * interval_ratio] = part
 
-            feeder_type, extra_slot = component_data.loc[part]['fdr'], 0
+            feeder_type, extra_slot = component_data.loc[part].fdr, 0
             extra_width = feeder_width[feeder_type][0] + feeder_width[feeder_type][1] - slot_interval
             while extra_width > 0:
                 extra_slot += 1
@@ -358,7 +403,7 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
                 extra_width -= slot_interval
 
             # 更新吸嘴信息
-            nozzle_pattern[idx] = component_data.loc[part]['nz']
+            nozzle_pattern[idx] = component_data.loc[part].nz
 
             # 更新头分配的先后顺序
             head_assign_indexes = np.array(best_assign_points).argsort().tolist()
@@ -376,11 +421,10 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
     if not optimal_nozzle_points:
         feeder_base, feeder_base_points = [-2] * max_slot_index, [0] * max_slot_index
         for _, feeder in feeder_data.iterrows():
-            slot, part = feeder['slot'], feeder['part']
-            part_index = component_data[component_data['part'] == part].index.tolist()[0]
+            part_index = component_data[component_data.part == feeder.part].index.tolist()[0]
 
             # 供料器基座分配位置和对应贴装点数
-            feeder_base[slot], feeder_base_points[slot] = part_index, feeder_division_points[part_index]
+            feeder_base[feeder.slot], feeder_base_points[feeder.slot] = part_index, feeder_division_points[part_index]
 
         # 前基座 TODO: 后基座
         for slot in range(max_slot_index // 2 - (max_head_index - 1) * interval_ratio):
@@ -397,12 +441,12 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
 
     # 更新供料器占位信息
     for _, data in feeder_data.iterrows():
-        feeder_base[data['slot']] = -1
+        feeder_base[data.slot] = -1
 
     for slot, feeder in enumerate(feeder_base):
         if feeder < 0:
             continue
-        part = component_data.loc[feeder]['part']
+        part = component_data.loc[feeder].part
 
         feeder_data.loc[len(feeder_data.index)] = [slot, part, 0]
 
@@ -414,20 +458,19 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
                      size=8)
 
         feeder_assign_range = []
-        for feeder in feeder_data.iterrows():
-            slot, part = feeder[1]['slot'], feeder[1]['part']
-            part_index = component_data[component_data['part'] == part].index.tolist()[0]
-            feeder_type = component_data.loc[part_index]['fdr']
+        for _, feeder in feeder_data.iterrows():
+            part_index = component_data[component_data.part == feeder.part].index.tolist()[0]
+            feeder_type = component_data.loc[part_index].fdr
             width = feeder_width[feeder_type][0] + feeder_width[feeder_type][1]
-            start = slotf1_pos[0] + slot_interval * (slot - 1) - slot_interval / 2
-            end = slotf1_pos[0] + slot_interval * (slot - 1) - slot_interval / 2 + width
+            start = slotf1_pos[0] + slot_interval * (feeder.slot - 1) - slot_interval / 2
+            end = slotf1_pos[0] + slot_interval * (feeder.slot - 1) - slot_interval / 2 + width
 
             rec_x = [start, end, end, start]
             rec_y = [slotf1_pos[1] - 40, slotf1_pos[1] - 40, slotf1_pos[1] + 10, slotf1_pos[1] + 10]
 
-            c = 'red' if feeder[1]['arg'] == 0 else 'black'        # 黑色表示已分配，红色表示新分配
-            plt.text(slotf1_pos[0] + slot_interval * (slot - 1), slotf1_pos[1] + 12,
-                     part + ': ' + str(feeder_points[part_index]), ha='center', size=7, rotation=90, color=c)
+            c = 'red' if feeder.arg == 0 else 'black'        # 黑色表示已分配，红色表示新分配
+            plt.text(slotf1_pos[0] + slot_interval * (feeder.slot - 1), slotf1_pos[1] + 12,
+                     feeder.part + ': ' + str(feeder_points[part_index]), ha='center', size=7, rotation=90, color=c)
 
             plt.fill(rec_x, rec_y, facecolor='yellow', alpha=0.4)
 
@@ -461,30 +504,28 @@ def feeder_allocate(component_data, pcb_data, feeder_data, figure=False):
 @timer_wrapper
 def feeder_base_scan(component_data, pcb_data, feeder_data, nozzle_pattern):
     feeder_assign_check = set()
-    for feeder in feeder_data.iterrows():
-        feeder_assign_check.add(feeder[1]['part'])
+    for _, feeder in feeder_data.iterrows():
+        feeder_assign_check.add(feeder.part)
 
     component_points = [0] * len(component_data)
-    for step in pcb_data.iterrows():
-        part = step[1]['part']
-        part_index = component_data[component_data['part'] == part].index.tolist()[0]
+    for i, data in pcb_data.iterrows():
+        part_index = component_data[component_data.part == data.part].index.tolist()[0]
 
         component_points[part_index] += 1
-        nozzle_type = component_data.loc[part_index]['nz']
+        nozzle_type = component_data.loc[part_index].nz
         if nozzle_type not in nozzle_limit.keys() or nozzle_limit[nozzle_type] <= 0:
             info = 'there is no available nozzle [' + nozzle_type + '] for the assembly process'
             raise ValueError(info)
 
     assert len(feeder_assign_check) == len(component_points) - component_points.count(0)    # 所有供料器均已分配槽位
     feeder_part = [-1] * max_slot_index
-    for feeder in feeder_data.iterrows():
-        part, slot = feeder[1]['part'], feeder[1]['slot']
-        part_index = component_data[component_data['part'] == part].index.tolist()
-        if len(part_index) != 1:
-            print('unregistered component: ', part, ' in slot', slot)
+    for _, data in feeder_data.iterrows():
+        component_index = component_data[component_data.part == data.part].index.tolist()
+        if len(component_index) != 1:
+            print('unregistered component: ', data.part, ' in slot', data.slot)
             continue
-        part_index = part_index[0]
-        feeder_part[slot] = part_index
+        component_index = component_index[0]
+        feeder_part[data.slot] = component_index
 
     component_result, cycle_result, feeder_slot_result = [], [], []  # 贴装点索引和拾取槽位优化结果
 
@@ -539,7 +580,7 @@ def feeder_base_scan(component_data, pcb_data, feeder_data, nozzle_pattern):
                             if scan_part[head] == -1 and part != -1 and component_points[part] > 0 and scan_part.count(
                                     part) < component_points[part]:
                                 # 2.匹配条件满足：不超过可用吸嘴数的限制
-                                nozzle = component_data.loc[part]['nz']
+                                nozzle = component_data.loc[part].nz
                                 if scan_nozzle_limit[nozzle] <= 0:
                                     continue
 
@@ -591,7 +632,7 @@ def feeder_base_scan(component_data, pcb_data, feeder_data, nozzle_pattern):
                         for head, nozzle in enumerate(nozzle_cycle):
                             if scan_part[head] == -1:
                                 continue
-                            if component_data.loc[scan_part[head]]['nz'] != nozzle and nozzle != '':
+                            if component_data.loc[scan_part[head]].nz != nozzle and nozzle != '':
                                 nozzle_counter += 2
 
                         # 下一周期（额外增加的吸嘴更换次数）
@@ -602,7 +643,7 @@ def feeder_base_scan(component_data, pcb_data, feeder_data, nozzle_pattern):
                                 prev_counter, new_counter = 0, 0
                                 if nozzle_cycle[head] != nozzle and nozzle_cycle[head] != '' and nozzle != '':
                                     prev_counter += 2
-                                if component_data.loc[scan_part[head]]['nz'] != nozzle and nozzle != '':
+                                if component_data.loc[scan_part[head]].nz != nozzle and nozzle != '':
                                     new_counter += 2
                                 nozzle_counter += new_counter - prev_counter
                         else:
@@ -612,7 +653,7 @@ def feeder_base_scan(component_data, pcb_data, feeder_data, nozzle_pattern):
                                 prev_counter, new_counter = 0, 0
                                 if nozzle_cycle[head] != nozzle and nozzle_cycle[head] != '' and nozzle != '':
                                     prev_counter += 2
-                                if component_data.loc[scan_part[head]]['nz'] != nozzle and nozzle != '':
+                                if component_data.loc[scan_part[head]].nz != nozzle and nozzle != '':
                                     new_counter += 2
                                 nozzle_counter += new_counter - prev_counter
 
@@ -708,7 +749,7 @@ def feeder_base_scan(component_data, pcb_data, feeder_data, nozzle_pattern):
             for head, component in enumerate(assigned_part):
                 if component == -1:
                     continue
-                cycle_nozzle[head] = component_data.loc[component]['nz']
+                cycle_nozzle[head] = component_data.loc[component].nz
 
             nozzle_mode.insert(nozzle_insert_cycle + 1, cycle_nozzle)
 

base_optimizer/optimizer_hybridgenetic.py

@@ -1,12 +1,4 @@
-import copy
-import random
-
-import matplotlib.pyplot as plt
-import numpy as np
-import pandas as pd
-
 from base_optimizer.optimizer_common import *
-from collections import defaultdict
 
 
 def dynamic_programming_cycle_path(cycle_placement, cycle_points):

base_optimizer/optimizer_interface.py

@@ -0,0 +1,96 @@
+# 用于提供对外接口
+from base_optimizer.optimizer_scanbased import *
+from base_optimizer.optimizer_celldivision import *
+from base_optimizer.optimizer_hybridgenetic import *
+from base_optimizer.optimizer_feederpriority import *
+from base_optimizer.optimizer_aggregation import *
+from base_optimizer.optimizer_twophase import *
+from base_optimizer.optimizer_mathmodel import *
+
+from base_optimizer.result_analysis import *
+
+
+def base_optimizer(machine_index, pcb_data, component_data, feeder_data=None, method='', hinter=False):
+
+    if method == 'cell_division':  # 基于元胞分裂的遗传算法
+        component_result, cycle_result, feeder_slot_result = optimizer_celldivision(pcb_data, component_data,
+                                                                                    hinter=False)
+        placement_result, head_sequence = greedy_placement_route_generation(component_data, pcb_data, component_result,
+                                                                            cycle_result, feeder_slot_result)
+    elif method == 'feeder_scan':  # 基于基座扫描的供料器优先算法
+        # 第1步：分配供料器位置
+        nozzle_pattern = feeder_allocate(component_data, pcb_data, feeder_data, figure=False)
+        # 第2步：扫描供料器基座，确定元件拾取的先后顺序
+        component_result, cycle_result, feeder_slot_result = feeder_base_scan(component_data, pcb_data, feeder_data,
+                                                                              nozzle_pattern)
+
+        # 第3步：贴装路径规划
+        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_for_route_generation(component_data, pcb_data, component_result,
+        #                                                                    cycle_result, feeder_slot_result)
+
+    elif method == 'hybrid_genetic':  # 基于拾取组的混合遗传算法
+        component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_hybrid_genetic(
+            pcb_data, component_data, hinter=False)
+
+    elif method == 'aggregation':  # 基于batch-level的整数规划 + 启发式算法
+        component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_aggregation(
+            component_data, pcb_data)
+    elif method == 'genetic_scanning':
+        component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_genetic_scanning(
+            component_data, pcb_data, hinter=False)
+    elif method == 'mip_model':
+        component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_mathmodel(
+            component_data, pcb_data, hinter=True)
+    elif method == "two_phase":
+        component_result, feeder_slot_result, cycle_result = gurobi_optimizer(pcb_data, component_data, feeder_data,
+                                                                              initial=True, partition=True,
+                                                                              reduction=True, hinter=hinter)
+
+        placement_result, head_sequence = scan_based_placement_route_generation(component_data, pcb_data,
+                                                                                component_result, cycle_result)
+    else:
+        raise 'method is not existed'
+
+    info = OptInfo()
+    assigned_nozzle = ['' if idx == -1 else component_data.loc[idx]['nz'] for idx in component_result[0]]
+    info.cycle_counter = sum(cycle_result)
+
+    for cycle in range(len(cycle_result)):
+        pick_slot = set()
+        for head in range(max_head_index):
+            idx = component_result[cycle][head]
+            if idx == -1:
+                continue
+
+            nozzle = component_data.loc[idx]['nz']
+            if nozzle != assigned_nozzle[head]:
+                if assigned_nozzle[head] != '':
+                    info.nozzle_change_counter += 1
+                assigned_nozzle[head] = nozzle
+
+            pick_slot.add(feeder_slot_result[cycle][head] - head * interval_ratio)
+        info.pickup_counter += len(pick_slot) * cycle_result[cycle]
+
+        pick_slot = list(pick_slot)
+        pick_slot.sort()
+        for idx in range(len(pick_slot) - 1):
+            info.pickup_movement += abs(pick_slot[idx + 1] - pick_slot[idx])
+
+    if hinter:
+        optimization_assign_result(component_data, pcb_data, component_result, cycle_result, feeder_slot_result,
+                                   nozzle_hinter=True, component_hinter=False, feeder_hinter=True)
+
+        print('----- Placement machine ' + str(machine_index) + ' ----- ')
+        print('-Cycle counter: {}'.format(info.cycle_counter))
+
+        print('-Nozzle change counter: {}'.format(info.nozzle_change_counter))
+        print('-Pick operation counter: {}'.format(info.pickup_counter))
+        print('-Pick movement: {}'.format(info.pickup_movement))
+        print('------------------------------ ')
+
+    # 估算贴装用时
+    info.placement_time = placement_time_estimate(component_data, pcb_data, component_result, cycle_result,
+                                                  feeder_slot_result, placement_result, head_sequence, hinter=False)
+    return info

base_optimizer/optimizer_mathmodel.py

@@ -0,0 +1,354 @@
+from base_optimizer.optimizer_common import *
+
+
+def list_range(start, end=None):
+    return list(range(start)) if end is None else list(range(start, end))
+
+
+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)
+
+    H = max_head_index
+    I = len(component_data)
+    S = len(component_data)
+    K = len(pcb_data)
+
+    nozzle_type, component_type = [], []
+    for _, data in component_data.iterrows():
+        if not data.nz in nozzle_type:
+            nozzle_type.append(data.nz)
+        component_type.append(data.part)
+
+    average_pos = 0
+    for _, data in pcb_data.iterrows():
+        average_pos += data.x
+    slot_start = int(round(average_pos / len(pcb_data) + stopper_pos[0] - slotf1_pos[0]) / slot_interval) + 1
+
+    r = 1
+    J = len(nozzle_type)
+    M = 10000
+    CompOfNozzle = [[0 for _ in range(J)] for _ in range(I)]    # Compatibility
+
+    component_point = [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
+
+    # objective related
+    g = mdl.addVars(list_range(K), vtype=GRB.BINARY)
+    d = mdl.addVars(list_range(K - 1), list_range(H), vtype=GRB.CONTINUOUS)
+    u = mdl.addVars(list_range(K), vtype=GRB.INTEGER)
+
+    d_plus = mdl.addVars(list_range(J), list_range(H), list_range(K - 1), vtype=GRB.CONTINUOUS)
+    d_minus = mdl.addVars(list_range(J), list_range(H), list_range(K - 1), vtype=GRB.CONTINUOUS)
+
+    e = mdl.addVars(list_range(-(H - 1) * r, S), list_range(K), vtype=GRB.BINARY)
+    f = mdl.addVars(list_range(S), list_range(I), vtype=GRB.BINARY, name='')
+    x = mdl.addVars(list_range(I), list_range(S), list_range(K), list_range(H), vtype=GRB.BINARY)
+    n = mdl.addVars(list_range(H), vtype=GRB.CONTINUOUS)
+
+    mdl.addConstrs(g[k] <= g[k + 1] for k in range(K - 1))
+
+    mdl.addConstrs(
+        quicksum(x[i, s, k, h] for i in range(I) for s in range(S)) <= g[k] for k in range(K) for h in range(H))
+
+    # nozzle no more than 1 for head h and cycle k
+    mdl.addConstrs(
+        quicksum(CompOfNozzle[i][j] * x[i, s, k, h] for i in range(I) for s in range(S) for j in range(J)) <= 1 for k in
+        range(K) for h in range(H))
+
+    # nozzle available number constraint
+    mdl.addConstrs(
+        quicksum(CompOfNozzle[i][j] * x[i, s, k, h] for i in range(I) for s in range(S) for h in range(H)) <= H for k in
+        range(K) for j in range(J))
+
+    # work completion
+    mdl.addConstrs(
+        quicksum(x[i, s, k, h] for s in range(S) for k in range(K) for h in range(H)) == component_point[i] for i in
+        range(I))
+
+    # nozzle change
+    mdl.addConstrs(quicksum(CompOfNozzle[i][j] * x[i, s, k, h] for i in range(I) for s in range(S)) - quicksum(
+        CompOfNozzle[i][j] * x[i, s, k + 1, h] for i in range(I) for s in range(S)) == d_plus[j, h, k] - d_minus[
+                       j, h, k] for k in range(K - 1) for j in range(J) for h in range(H))
+
+    mdl.addConstrs(
+        2 * d[k, h] == quicksum(d_plus[j, h, k] for j in range(J)) + quicksum(d_minus[j, h, k] for j in range(J)) for k
+        in range(K - 1) for h in range(H))
+
+    mdl.addConstrs(n[h] == quicksum(d[k, h] for k in range(K - 1)) - 0.5 for h in range(H))
+
+    # simultaneous pick
+    for s in range(-(H - 1) * r, S):
+        rng = list(range(max(0, -math.floor(s / r)), min(H, math.ceil((S - s) / r))))
+        for k in range(K):
+            mdl.addConstr(quicksum(x[i, s + h * r, k, h] for h in rng for i in range(I)) <= M * e[s, k], name='')
+            mdl.addConstr(quicksum(x[i, s + h * r, k, h] for h in rng for i in range(I)) >= e[s, k], name='')
+    # pickup movement
+    mdl.addConstrs(
+        u[k] >= s1 * e[s1, k] - s2 * e[s2, k] for s1 in range(-(H - 1) * r, S) for s2 in 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 i in range(I)) <= 1 for s in range(S))
+    mdl.addConstrs(
+        quicksum(x[i, 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))
+    mdl.addConstrs(
+        quicksum(x[i, s, k, h] for h in range(H) for k in range(K)) <= M * f[s, i] for i in range(I) for s in
+        range(S))
+
+    # objective
+    t_c, t_n, t_p, t_m = 2, 6, 1, 0.1
+    mdl.setObjective(t_c * quicksum(g[k] for k in range(K)) + t_n * quicksum(
+        d[k, h] for h in range(H) for k in range(K - 1)) + t_p * quicksum(
+        e[s, k] for s in range(-(H - 1) * r, S) for k in range(K)) + t_m * quicksum(u[k] for k in range(K)),
+                     GRB.MINIMIZE)
+
+    mdl.optimize()
+
+    component_result, cycle_result, feeder_slot_result = [], [], []
+    for k in range(K):
+        if abs(g[k].x) < 1e-6:
+            continue
+
+        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):
+                for s in range(S):
+                    if abs(x[i, s, k, h].x) > 1e-6:
+                        component_result[-1][h] = i
+                        feeder_slot_result[-1][h] = slot_start + s * interval_ratio - 1
+    if hinter:
+        print(component_result)
+        print(feeder_slot_result)
+
+    return component_result, cycle_result, feeder_slot_result
+
+
+def place_route_model(component_data, pcb_data, component_result, feeder_slot_result, figure=False, hinter=True):
+    mdl = Model('place_route')
+    mdl.setParam('Seed', 0)
+    mdl.setParam('OutputFlag', hinter)  # set whether output the debug information
+    # mdl.setParam('TimeLimit', 20)
+
+    component_type = []
+    for _, data in component_data.iterrows():
+        component_type.append(data.part)
+
+    pos = []
+    for _, data in pcb_data.iterrows():
+        pos.append([data.x + stopper_pos[0], data.y + stopper_pos[1]])
+
+    I, P, H = len(component_data), len(pcb_data), max_head_index
+    A = []
+    for h1 in range(H):
+        for h2 in range(H):
+            if h1 == h2:
+                continue
+            A.append([h1, h2])
+    K = len(component_result)
+
+    CompOfPoint = [[0 for _ in range(P)] for _ in range(I)]
+    for row, data in pcb_data.iterrows():
+        idx = component_type.index(data.part)
+        CompOfPoint[idx][row] = 1
+
+    d_FW, d_PL, d_BW = np.zeros([P, K, H]), np.zeros([P, P, len(A)]), np.zeros([P, K, H])
+    for k in range(K):
+        min_slot, max_slot = float('inf'), float('-inf')
+        for h in range(H):
+            if feeder_slot_result[k][h] == -1:
+                continue
+            min_slot = min(min_slot, feeder_slot_result[k][h] - h * interval_ratio)
+            max_slot = max(max_slot, feeder_slot_result[k][h] - h * interval_ratio)
+
+        for p in range(P):
+            for h in range(H):
+                d_FW[p, k, h] = max(
+                    abs(slotf1_pos[0] + (max_slot - 1) * slot_interval - pos[p][0] + h * head_interval),
+                    abs(slotf1_pos[1] - pos[p][1]))
+
+                d_BW[p, k, h] = max(
+                    abs(slotf1_pos[0] + (min_slot - 1) * slot_interval - pos[p][0] + h * head_interval),
+                    abs(slotf1_pos[1] - pos[p][1]))
+
+    for p in range(P):
+        for q in range(P):
+            for idx, arc in enumerate(A):
+                h1, h2 = arc
+                d_PL[p, q, idx] = max(abs(pos[p][0] - pos[q][0] - (h1 - h2) * head_interval), abs(pos[p][1] - pos[q][1]))
+
+    w = mdl.addVars(list_range(P), list_range(P), list_range(K), list_range(len(A)), vtype=GRB.BINARY)
+    y = mdl.addVars(list_range(P), list_range(K), list_range(H), vtype=GRB.BINARY)
+    z = mdl.addVars(list_range(P), list_range(K), list_range(H), vtype=GRB.BINARY)
+
+    def A_from(h):
+        res = []
+        for idx, arc in enumerate(A):
+            if arc[0] == h:
+                res.append(idx)
+        return res
+
+    def A_to(h):
+        res = []
+        for idx, arc in enumerate(A):
+            if arc[1] == h:
+                res.append(idx)
+        return res
+
+    def A_contain(h):
+        res = []
+        for idx, arc in enumerate(A):
+            if h in arc:
+                res.append(idx)
+        return res
+
+    # constraints on component assignment type, assigned points cannot conflict with the corresponding component type
+    for k in range(K):
+        for h in range(H):
+            if component_result[k][h] == -1:
+                # no components on the head
+                mdl.addConstr(quicksum(w[p, q, k, a] for a in A_contain(h) for q in range(P) for p in range(P)) == 0)
+            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))) / 2 <= CompOfPoint[component_result[k][h]][p] for
+                               p in range(P))
+
+    # each head corresponds to a maximum of one point in each cycle
+    mdl.addConstrs(
+        quicksum(w[p, q, k, a] for p in range(P) for q in range(P) for a in A_contain(h)) <= 2 for k in range(K) for h
+        in range(H))
+
+    mdl.addConstrs(
+        quicksum((y[p, k, h] + z[p, k, h]) for p in range(P)) <= 1 for k in range(K) for h in
+        range(H))
+
+    # task continuity (for the same point the entering head and the leaving head should be same)
+    mdl.addConstrs(quicksum(w[p, q, k, a] for p in range(P) for a in A_to(h)) + y[q, k, h] == quicksum(
+        w[q, p, k, a] for p in range(P) for a in A_from(h)) + z[q, k, h] for k in range(K) for h in range(H) for q in
+                   range(P))
+
+    mdl.addConstrs(
+        y[p, k, h] <= quicksum(w[p, q, k, a] for q in range(P) for a in A_from(h)) for h in range(H) for p in
+        range(P) for k in range(K))
+
+    mdl.addConstrs(
+        z[p, k, h] <= quicksum(w[q, p, k, a] for q in range(P) for a in A_to(h)) for h in range(H) for p in
+        range(P) for k in range(K))
+
+    # one arrival point per cycle
+    mdl.addConstrs(quicksum(y[p, k, h] for p in range(P) for h in range(H)) == 1 for k in range(K))
+
+    # one departure point per cycle
+    mdl.addConstrs(quicksum(z[p, k, h] for p in range(P) for h in range(H)) == 1 for k in range(K))
+
+    # one enter edge per point
+    mdl.addConstrs(quicksum(y[q, k, h] for h in range(H) for k in range(K)) + quicksum(
+        w[p, q, k, a] for p in range(P) for a in range(len(A)) for k in range(K)) == 1 for q in range(P))
+
+    # one leaving edge per point
+    mdl.addConstrs(quicksum(z[q, k, h] for h in range(H) for k in range(K)) + quicksum(
+        w[q, p, k, a] for p in range(P) for a in range(len(A)) for k in range(K)) == 1 for q in range(P))
+
+    # subtour eliminate constraint
+    n = mdl.addVars(list_range(P), vtype=GRB.CONTINUOUS)
+    m = mdl.addVars(list_range(P), vtype=GRB.CONTINUOUS)
+    v = mdl.addVars(list_range(P), list_range(P), vtype=GRB.CONTINUOUS)
+
+    mdl.addConstrs(
+        m[p] + quicksum(v[p, q] for q in range(P)) - n[p] - quicksum(v[q, p] for q in range(P)) == 1 for p in range(P))
+
+    mdl.addConstrs(
+        v[p, q] <= (P - K + 1) * quicksum(w[p, q, k, a] for a in range(len(A)) for k in range(K)) for p in range(P) for
+        q in range(P))
+
+    mdl.addConstrs(n[p] <= (P - K + 1) * quicksum(y[p, k, h] for h in range(H) for k in range(K)) for p in range(P))
+    mdl.addConstrs(m[p] <= (P - K + 1) * quicksum(z[p, k, h] for h in range(H) for k in range(K)) for p in range(P))
+
+    # objective
+    mdl.setObjective(
+        quicksum(d_FW[p, k, h] * y[p, k, h] for p in range(P) for k in range(K) for h in range(H)) + quicksum(
+            d_PL[p, q, a] * w[p, q, k, a] for k in range(K) for p in range(P) for q in range(P) for a in
+            range(len(A))) + quicksum(d_BW[p, k, h] * z[p, k, h] for p in range(P) for k in range(K) for h in range(H)),
+        GRB.MINIMIZE)
+
+    mdl.optimize()
+    if figure:
+        for k in range(K):
+            plt.scatter([p[0] for p in pos[0:8]], [p[1] for p in pos[0:8]], color='red')
+            plt.scatter([p[0] for p in pos[8:]], [p[1] for p in pos[8:]], color='blue')
+            for p in range(P):
+                for q in range(P):
+                    for idx, arc in enumerate(A):
+                        if abs(w[p, q, k, idx].x) > 1e-6:
+                            h1, h2 = arc
+                            plt.plot([pos[p][0] - h1 * head_interval, pos[q][0] - h2 * head_interval],
+                                     [pos[p][1], pos[q][1]], linestyle='-.', color='black', linewidth=1)
+                            plt.text(pos[p][0] - h1 * head_interval, pos[p][1], 'H%d' % (h1 + 1), ha='center',
+                                     va='bottom', size=10)
+
+                for h in range(H):
+                    if abs(y[p, k, h].x) > 1e-6:
+                        plt.plot([pos[p][0] - h * head_interval, 500], [pos[p][1], 100], linestyle='-.', color='black',
+                                 linewidth=1)
+                        plt.text(pos[p][0] - h * head_interval, pos[p][1], 'H%d' % (h + 1), ha='center', va='bottom',
+                                 size=10)
+
+                for h in range(H):
+                    if abs(z[p, k, h].x) > 1e-6:
+                        plt.plot([pos[p][0] - h * head_interval, 900], [pos[p][1], 100], linestyle='-.', color='black',
+                                 linewidth=1)
+                        plt.text(pos[p][0] - h * head_interval, pos[p][1], 'H%d' % (h + 1), ha='center', va='bottom',
+                                 size=10)
+            plt.show()
+
+    # convert model result into standard form
+    placement_result, head_sequence = [[-1 for _ in range(H)] for _ in range(K)], [[] for _ in
+                                                                                   range(K)]
+    for k in range(K):
+        arc_list = []
+        for p in range(P):
+            for q in range(P):
+                for idx, arc in enumerate(A):
+                    if abs(w[p, q, k, idx].x) > 1e-6:
+                        plt.plot([pos[p][0], pos[q][0]], [pos[p][1], pos[q][1]], linestyle='-.', color='black',
+                                 linewidth=1)
+                        placement_result[k][arc[0]], placement_result[k][arc[1]] = p, q
+                        arc_list.append(arc)
+
+        head, idx = -1, 0
+        for p in range(P):
+            for h in range(H):
+                if abs(y[p, k, h].x) > 1e-6:
+                    head = h
+
+        while idx < len(arc_list):
+            for i, arc in enumerate(arc_list):
+                if arc[0] == head:
+                    head_sequence[k].append(head)
+                    head = arc[1]
+                    idx += 1
+                    break
+        head_sequence[k].append(head)
+
+    return placement_result, head_sequence
+
+
+@timer_wrapper
+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)
+    return component_result, cycle_result, feeder_slot_result, placement_result, head_sequence

base_optimizer/optimizer_scanbased.py

@@ -1,4 +1,3 @@
-import itertools
 from base_optimizer.optimizer_common import *
 
 

base_optimizer/optimizer_twophase.py

@@ -0,0 +1,855 @@
+from base_optimizer.optimizer_common import *
+
+
+def list_range(start, end=None):
+    return list(range(start)) if end is None else list(range(start, end))
+
+
+@timer_wrapper
+def gurobi_optimizer(pcb_data, component_data, feeder_data, reduction=True, partition=True, initial=False, hinter=True):
+    # data preparation: convert data to index
+    component_list, nozzle_list = defaultdict(int), defaultdict(int)
+    cpidx_2_part, nzidx_2_nozzle, cpidx_2_nzidx = {}, {}, {}
+    arg_slot_rng = None if len(feeder_data) == 0 else [feeder_data.iloc[0].slot, feeder_data.iloc[-1].slot]
+    for idx, data in component_data.iterrows():
+        part, nozzle = data.part, data.nz
+
+        cpidx_2_part[idx] = part
+        nz_key = [key for key, val in nzidx_2_nozzle.items() if val == nozzle]
+
+        nz_idx = len(nzidx_2_nozzle) if len(nz_key) == 0 else nz_key[0]
+        nzidx_2_nozzle[nz_idx] = nozzle
+
+        component_list[part] = 0
+        cpidx_2_nzidx[idx] = nz_idx
+
+    for _, data in pcb_data.iterrows():
+        idx = component_data[component_data.part == data.part].index.tolist()[0]
+        nozzle = component_data.loc[idx].nz
+
+        nozzle_list[nozzle] += 1
+        component_list[data.part] += 1
+
+    part_feederbase = defaultdict(int)
+    if feeder_data is not None:
+        for _, data in feeder_data.iterrows():
+            idx = -1
+            for idx, part_ in cpidx_2_part.items():
+                if data.part == part_:
+                    break
+            assert idx != -1
+            part_feederbase[idx] = data.slot  # part index - slot
+
+    if not reduction:
+        ratio = 2  # 直接导入飞达数据时，采用正常吸杆间隔
+    else:
+        if len(component_list) <= 1.5 * max_head_index:
+            ratio = 1
+        else:
+            ratio = 2
+    I, J = len(cpidx_2_part.keys()), len(nzidx_2_nozzle.keys())
+    # === determine the hyper-parameter of L ===
+    # first phase: calculate the number of heads for each type of nozzle
+    nozzle_heads = defaultdict(int)
+    for nozzle in nozzle_list.keys():
+        nozzle_heads[nozzle] = 1
+
+    while sum(nozzle_heads.values()) != max_head_index:
+        max_cycle_nozzle = None
+
+        for nozzle, head_num in nozzle_heads.items():
+            if max_cycle_nozzle is None or nozzle_list[nozzle] / head_num > nozzle_list[max_cycle_nozzle] / \
+                    nozzle_heads[max_cycle_nozzle]:
+                max_cycle_nozzle = nozzle
+
+        assert max_cycle_nozzle is not None
+        nozzle_heads[max_cycle_nozzle] += 1
+
+    nozzle_comp_points = defaultdict(list)
+    for part, points in component_list.items():
+        idx = component_data[component_data.part == part].index.tolist()[0]
+        nozzle = component_data.loc[idx].nz
+        nozzle_comp_points[nozzle].append([part, points])
+
+    level = 1 if len(component_list) == 1 or len(component_list) % max_head_index == 0 else 2
+    part_assignment, cycle_assignment = [], []
+
+    def aux_func(info):
+        return max(map(lambda points: max([p[1] for p in points]), info))
+
+    pre_objbst, pre_changetime = None, None
+    def terminate_condition(mdl, where):
+        if where == GRB.Callback.MIP:
+            objbst = mdl.cbGet(GRB.Callback.MIP_OBJBST)
+            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 > 45:
+                    # pass
+                    mdl.terminate()
+            else:
+                pre_changetime = changetime
+
+            pre_objbst = objbst
+
+    def recursive_assign(assign_points, nozzle_compo_points, cur_level, total_level) -> int:
+        def func(points):
+            return map(lambda points: max([p[1] for p in points]), points)
+
+        if cur_level > total_level and sum(func(nozzle_compo_points.values())) == 0:
+            return 0
+        elif assign_points <= 0 and cur_level == 1:
+            return -1  # backtrack
+        elif assign_points <= 0 or cur_level > total_level:
+            return 1  # fail
+
+        nozzle_compo_points_cpy = copy.deepcopy(nozzle_compo_points)
+        prev_assign = 0
+        for part in part_assignment[cur_level - 1]:
+            if part != -1:
+                prev_assign += 1
+
+        head_idx = 0
+        for nozzle, head in nozzle_heads.items():
+            while head:
+                min_idx = -1
+                for idx, (part, points) in enumerate(nozzle_compo_points_cpy[nozzle]):
+                    if points >= assign_points and (
+                            min_idx == -1 or points < nozzle_compo_points_cpy[nozzle][min_idx][1]):
+                        min_idx = idx
+                part_assignment[cur_level - 1][head_idx] = -1 if min_idx == -1 else \
+                    nozzle_compo_points_cpy[nozzle][min_idx][0]
+                if min_idx != -1:
+                    nozzle_compo_points_cpy[nozzle][min_idx][1] -= assign_points
+                head -= 1
+                head_idx += 1
+
+        cycle_assignment[cur_level - 1] = assign_points
+        for part in part_assignment[cur_level - 1]:
+            if part != -1:
+                prev_assign -= 1
+
+        if prev_assign == 0:
+            res = 1
+        else:
+            points = min(len(pcb_data) // max_head_index + 1, aux_func(nozzle_compo_points_cpy.values()))
+            res = recursive_assign(points, nozzle_compo_points_cpy, cur_level + 1, total_level)
+        if res == 0:
+            return 0
+        elif res == 1:
+            # All cycles have been completed, but there are still points left to be allocated
+            return recursive_assign(assign_points - 1, nozzle_compo_points, cur_level, total_level)
+
+    # second phase: (greedy) recursive search to assign points for each cycle set and obtain an initial solution
+    while True:
+        part_assignment = [[-1 for _ in range(max_head_index)] for _ in range(level)]
+        cycle_assignment = [-1 for _ in range(level)]
+        points = min(len(pcb_data) // max_head_index + 1, max(component_list.values()))
+        if recursive_assign(points, nozzle_comp_points, 1, level) == 0:
+            break
+        level += 1
+
+    weight_cycle, weight_nz_change, weight_pick = 2, 3, 2
+
+    L = len(cycle_assignment) if partition else len(pcb_data)
+    S = ratio * I 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):
+        for j in range(J):
+            HC[i][j] = 1 if cpidx_2_nzidx[i] == j else 0
+
+    mdl = Model('SMT')
+    mdl.setParam('Seed', 0)
+    mdl.setParam('OutputFlag', hinter)  # set whether output the debug information
+    mdl.setParam('TimeLimit', 3600 * 3)
+    mdl.setParam('PoolSearchMode', 2)
+    mdl.setParam('PoolSolutions', 3e2)
+    mdl.setParam('PoolGap', 1e-4)
+    # mdl.setParam('MIPFocus', 2)
+    # 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)
+
+    # === Decision Variables ===
+    x = mdl.addVars(list_range(I), list_range(S), list_range(max_head_index), list_range(L), vtype=GRB.BINARY, name='x')
+    y = mdl.addVars(list_range(I), list_range(max_head_index), list_range(L), vtype=GRB.BINARY, name='y')
+    v = mdl.addVars(list_range(S), list_range(max_head_index), list_range(L), vtype=GRB.BINARY, name='v')
+
+    c = mdl.addVars(list_range(I), list_range(max_head_index), list_range(L), vtype=GRB.INTEGER, name='c')
+
+    mdl.addConstrs(
+        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():
+            for s in range(S):
+                f[s, i] = mdl.addVar(vtype=GRB.BINARY, name='f_' + str(s) + '_' + str(i))
+        else:
+            for s in range(S):
+                f[s, i] = 1 if part_feederbase[i] == s + arg_slot_rng[0] else 0
+
+    p = mdl.addVars(list_range(-(max_head_index - 1) * ratio, S), list_range(L), vtype=GRB.BINARY, name='p')
+    z = mdl.addVars(list_range(J), list_range(max_head_index), list_range(L), vtype=GRB.BINARY)
+
+    d = mdl.addVars(list_range(L), list_range(max_head_index), vtype=GRB.INTEGER, name='d')
+    d_plus = mdl.addVars(list_range(J), list_range(max_head_index), list_range(L), vtype=GRB.INTEGER,
+                         name='d_plus')
+    d_minus = mdl.addVars(list_range(J), list_range(max_head_index), list_range(L), vtype=GRB.INTEGER,
+                          name='d_minus')
+
+    max_cycle = math.ceil(len(pcb_data) / max_head_index)
+    PU = mdl.addVars(list_range(-(max_head_index - 1) * ratio, S), list_range(L), vtype=GRB.INTEGER, name='PU')
+    WL = mdl.addVars(list_range(L), vtype=GRB.INTEGER, ub=max_cycle, name='WL')
+    NC = mdl.addVars(list_range(max_head_index), vtype=GRB.INTEGER, name='NC')
+
+    part_2_cpidx = defaultdict(int)
+    for idx, part in cpidx_2_part.items():
+        part_2_cpidx[part] = idx
+
+    if initial:
+        # initial some variables to speed up the search process
+        # ensure the priority of the workload assignment
+        cycle_index = sorted(range(len(cycle_assignment)), key=lambda k: cycle_assignment[k], reverse=True)
+        part_list = []
+
+        for cycle in cycle_index:
+            cycle_part = part_assignment[cycle]
+            for part in cycle_part:
+                if part != -1 and part not in part_list:
+                    part_list.append(part)
+        slot = 0
+        for part in part_list:
+            if feeder_data is not None:
+                while slot in feeder_data.keys():
+                    slot += 1  # skip assigned feeder slot
+
+            if part_2_cpidx[part] in part_feederbase.keys():
+                continue
+
+            part_feederbase[part_2_cpidx[part]] = slot
+            # f[slot, part_2_cpidx[part]].Start = 1
+            slot += 1
+
+        for idx, cycle in enumerate(cycle_index):
+            WL[idx].Start = cycle_assignment[cycle]
+            for h in range(max_head_index):
+                part = part_assignment[cycle][h]
+                if part == -1:
+                    continue
+                i = part_2_cpidx[part]
+                y[i, h, idx].Start = 1
+                v[part_feederbase[i], h, idx].Start = 1
+
+    # === Objective ===
+    mdl.setObjective(weight_cycle * quicksum(WL[l] for l in range(L)) + weight_nz_change * quicksum(
+        NC[h] for h in range(max_head_index)) + weight_pick * quicksum(
+        PU[s, l] for s in range(-(max_head_index - 1) * ratio, S) for l in range(L)))
+
+    # === Constraint ===
+    if not partition:
+        mdl.addConstrs(WL[l] <= 1 for l in range(L))
+
+    # work completion
+    # mdl.addConstrs(c[i, h, l] == WL[l] * y[i, h, l] for i in range(I) for h in range(max_head_index) for l in range(L))
+    mdl.addConstrs(
+        c[i, h, l] <= max_cycle * y[i, h, l] for i in range(I) for h in range(max_head_index) for l in range(L))
+    mdl.addConstrs(c[i, h, l] <= WL[l] for i in range(I) for h in range(max_head_index) for l in range(L))
+    mdl.addConstrs(
+        c[i, h, l] >= WL[l] - max_cycle * (1 - y[i, h, l]) for i in range(I) for h in range(max_head_index) for l in
+        range(L))
+
+    mdl.addConstrs(
+        quicksum(c[i, h, l] for h in range(max_head_index) for l in range(L)) == component_list[cpidx_2_part[i]] for i
+        in range(I))
+
+    # variable constraint
+    mdl.addConstrs(quicksum(y[i, h, l] for i in range(I)) <= 1 for h in range(max_head_index) for l in range(L))
+
+    # simultaneous pick
+    for s in range(-(max_head_index - 1) * ratio, S):
+        rng = list(range(max(0, -math.floor(s / ratio)), min(max_head_index, math.ceil((S - s) / ratio))))
+        for l in range(L):
+            mdl.addConstr(quicksum(v[s + h * ratio, h, l] for h in rng) <= max_head_index * p[s, l])
+            mdl.addConstr(quicksum(v[s + h * ratio, h, l] for h in rng) >= p[s, l])
+
+    # mdl.addConstrs(PU[s, l] == p[s, l] * WL[l] for s in range(-(max_head_index - 1) * ratio, S) for l in range(L))
+    mdl.addConstrs(PU[s, l] <= max_cycle * p[s, l] for s in range(-(max_head_index - 1) * ratio, S) for l in range(L))
+    mdl.addConstrs(PU[s, l] <= WL[l] for s in range(-(max_head_index - 1) * ratio, S) for l in range(L))
+    mdl.addConstrs(
+        PU[s, l] >= WL[l] - max_cycle * (1 - p[s, l]) for s in range(-(max_head_index - 1) * ratio, S) for l in
+        range(L))
+
+    # nozzle change
+    mdl.addConstrs(
+        z[j, h, l] - z[j, h, l + 1] == d_plus[j, h, l] - d_minus[j, h, l] for l in range(L - 1) for j in range(J) for h
+        in range(max_head_index))
+
+    mdl.addConstrs(z[j, h, 0] - z[j, h, L - 1] == d_plus[j, h, L - 1] - d_minus[j, h, L - 1] for j in range(J) for h
+        in range(max_head_index))
+
+    mdl.addConstrs(
+        2 * d[l, h] == quicksum(d_plus[j, h, l] for j in range(J)) + quicksum(d_minus[j, h, l] for j in range(J)) for l
+        in range(L - 1) for h in range(max_head_index))
+
+    mdl.addConstrs(2 * d[L - 1, h] == quicksum(d_plus[j, h, L - 1] for j in range(J)) + quicksum(
+        d_minus[j, h, L - 1] for j in range(J)) for h in range(max_head_index))
+
+    mdl.addConstrs(NC[h] == quicksum(d[l, h] for l in range(L)) for h in range(max_head_index))
+
+    mdl.addConstrs(quicksum(y[i, h, l] for i in range(I) for h in range(max_head_index)) * M >= WL[l] for l in range(L))
+
+    # nozzle-component compatibility
+    mdl.addConstrs(
+        y[i, h, l] <= quicksum(HC[i][j] * z[j, h, l] for j in range(J)) for i in range(I) for h in range(max_head_index)
+        for l in range(L))
+
+    # available number of feeder
+    mdl.addConstrs(quicksum(f[s, i] for s in range(S)) <= 1 for i in range(I))
+
+    # available number of nozzle
+    mdl.addConstrs(quicksum(z[j, h, l] for h in range(max_head_index)) <= max_head_index for j in range(J) for l in range(L))
+
+    # upper limit for occupation for feeder slot
+    mdl.addConstrs(quicksum(f[s, i] for i in range(I)) <= 1 for s in range(S))
+    mdl.addConstrs(
+        quicksum(v[s, h, l] for s in range(S)) >= quicksum(y[i, h, l] for i in range(I)) for h in range(max_head_index)
+        for l in range(L))
+
+    # others
+    mdl.addConstrs(quicksum(z[j, h, l] for j in range(J)) <= 1 for h in range(max_head_index) for l in range(L))
+    # mdl.addConstrs(
+    #     quicksum(x[i, s, h, l] for h in range(max_head_index) for l in range(L)) >= f[s, i] for i in range(I)
+    #     for s in range(S))
+    # mdl.addConstrs(
+    #     quicksum(x[i, s, h, l] for h in range(max_head_index) for l in range(L)) <= M * f[s, i] for i in
+    #     range(I) for s in range(S))
+
+    mdl.addConstrs(
+        f[s, i] >= x[i, s, h, l] for s in range(S) for i in range(I) for h in range(max_head_index) for l in range(L))
+
+    mdl.addConstrs(
+        quicksum(x[i, s, h, l] for h in range(max_head_index) for l in range(L)) >= f[s, i] for s in
+        range(S) for i in range(I))
+
+    # the constraints to speed up the search process
+    mdl.addConstrs(
+        quicksum(x[i, s, h, l] for i in range(I) for s in range(S)) <= 1 for h in range(max_head_index) for l
+        in range(L))
+
+    if reduction:
+        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(WL[l] >= WL[l + 1] for l in range(L - 1))
+        # 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))
+
+    mdl.addConstrs(
+        x[i, s, h, l] >= y[i, h, l] + v[s, h, l] - 1 for i in range(I) for s in range(S) for h in range(max_head_index)
+        for l in range(L))
+    mdl.addConstrs(
+        x[i, s, h, l] <= y[i, h, l] for i in range(I) for s in range(S) for h in range(max_head_index)
+        for l in range(L))
+
+    mdl.addConstrs(
+        x[i, s, h, l] <= v[s, h, l] for i in range(I) for s in range(S) for h in range(max_head_index)
+        for l in range(L))
+
+    # === search process ===
+    mdl.update()
+    # mdl.write('mdl.lp')
+    if hinter:
+        print('num of constrs: ', str(len(mdl.getConstrs())), ', num of vars: ', str(len(mdl.getVars())))
+
+    mdl.optimize(terminate_condition)
+
+    # === result generation ===
+    nozzle_assign, component_assign = [], []
+    feeder_assign, cycle_assign = [], []
+    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)
+        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])
+
+        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])]
+
+        min_dist, solution_number = None, -1
+        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)
+            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
+
+            for h in range(max_head_index):
+                for l in range(L):
+                    if abs(WL[l].Xn) <= 1e-4:
+                        continue
+
+                    pos_list = []
+                    for i in range(I):
+                        if abs(y[i, h, l].Xn) <= 1e-4:
+                            continue
+
+                        for _ in range(round(WL[l].Xn)):
+                            pos_list.append(component_pos[i][pos_counter[i]])
+                            pos_counter[i] += 1
+
+                        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
+
+            if min_dist is None or dist < min_dist:
+                min_dist = dist
+                solution_number = sol_counter
+
+        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)])
+
+            cycle_assign.append(round(WL[l].Xn))
+            if abs(WL[l].Xn) <= 1e-4:
+                continue
+
+            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 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
+            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 l in range(L):
+            if abs(WL[l].Xn) <= 1e-4:
+                continue
+
+            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)
+
+        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
+
+
+def scan_based_placement_route_generation(component_data, pcb_data, component_assign, cycle_assign):
+    placement_result, head_sequence_result = [], []
+
+    mount_point_pos, mount_point_index, mount_point_angle, mount_point_part = [], [], [], []
+    for i, data in pcb_data.iterrows():
+        component_index = component_data[component_data.part == data.part].index.tolist()[0]
+        # 记录贴装点序号索引和对应的位置坐标
+        mount_point_index.append(i)
+        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)
+
+    ref_pos_y = min(mount_point_pos, key=lambda x: x[1])[1]
+    for cycle_index, component_cycle in enumerate(component_assign):
+        for _ in range(cycle_assign[cycle_index]):
+            min_dist = None
+            tmp_assigned_placement, tmp_assigned_head_seq = [], []
+            tmp_mount_point_pos, tmp_mount_point_index = [], []
+            for search_dir in range(3):     # 不同的搜索方向，贴装头和起始点的选取方法各不相同
+                if search_dir == 0:
+                    # 从左向右搜索
+                    searchPoints = np.arange(lBoundary, (lBoundary + rBoundary) / 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)
+                    head_range = list(range(max_head_index - 1, -1, -1))
+                else:
+                    # 从中间向两边搜索
+                    searchPoints = np.arange(lBoundary, rBoundary, search_step / 2)
+                    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 startPoint in searchPoints:
+                    mount_point_pos_cpy, mount_point_index_cpy = copy.deepcopy(mount_point_pos), copy.deepcopy(
+                        mount_point_index)
+                    mount_point_angle_cpy = copy.deepcopy(mount_point_angle)
+
+                    assigned_placement = [-1] * max_head_index
+                    assigned_mount_point = [[0, 0]] * max_head_index
+                    assigned_mount_angle = [0] * max_head_index
+                    head_counter, point_index = 0, -1
+                    for head_index in head_range:
+                        if head_counter == 0:
+                            component_index = component_assign[cycle_index][head_index]
+
+                            if component_index == -1:
+                                continue
+
+                            min_horizontal_distance = None
+                            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(
+                                    mount_point_pos_cpy[index][1] - ref_pos_y)
+
+                                if min_horizontal_distance is None or horizontal_distance < min_horizontal_distance:
+                                    min_horizontal_distance = horizontal_distance
+                                    point_index = index
+                        else:
+                            point_index = -1
+                            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
+                            for index, mount_index in enumerate(mount_point_index_cpy):
+                                if mount_point_part[mount_index] != next_comp_index:
+                                    continue
+
+                                point_pos = [[mount_point_pos_cpy[index][0] - head_index * head_interval,
+                                              mount_point_pos_cpy[index][1]]]
+
+                                cheby_distance, euler_distance = 0, 0
+                                for next_head in range(max_head_index):
+                                    if assigned_placement[next_head] == -1:
+                                        continue
+                                    point_pos.append(assigned_mount_point[next_head].copy())
+                                    point_pos[-1][0] -= next_head * head_interval
+
+                                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)
+
+                                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
+
+                        if point_index == -1:
+                            continue
+
+                        head_counter += 1
+
+                        assigned_placement[head_index] = mount_point_index_cpy[point_index]
+                        assigned_mount_point[head_index] = mount_point_pos_cpy[point_index].copy()
+                        assigned_mount_angle[head_index] = mount_point_angle_cpy[point_index]
+
+                        mount_point_index_cpy.pop(point_index)
+                        mount_point_pos_cpy.pop(point_index)
+                        mount_point_angle_cpy.pop(point_index)
+
+                    dist, head_seq = dynamic_programming_cycle_path(assigned_placement, assigned_mount_point,
+                                                                    assigned_mount_angle)
+
+                    if min_dist is None or dist < min_dist:
+                        tmp_mount_point_pos, tmp_mount_point_index = mount_point_pos_cpy, mount_point_index_cpy
+                        tmp_assigned_placement, tmp_assigned_head_seq = assigned_placement, head_seq
+                        min_dist = dist
+
+            mount_point_pos, mount_point_index = tmp_mount_point_pos, tmp_mount_point_index
+
+            placement_result.append(tmp_assigned_placement)
+            head_sequence_result.append(tmp_assigned_head_seq)
+
+    return placement_result, head_sequence_result
+    # return placement_route_relink_heuristic(component_data, pcb_data, placement_result, head_sequence_result)
+
+
+def placement_route_relink_heuristic(component_data, pcb_data, placement_result, head_sequence_result, hinter=True):
+    mount_point_pos, mount_point_angle, mount_point_index, mount_point_part = [], [], [], []
+    for i, data in pcb_data.iterrows():
+        component_index = component_data[component_data.part == data.part].index.tolist()[0]
+        # 记录贴装点序号索引和对应的位置坐标
+        mount_point_index.append(i)
+        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)
+
+    cycle_length, cycle_average_pos = [], []
+    for cycle, placement in enumerate(placement_result):
+        prev_pos, prev_angle = None, None
+        cycle_pos_list = []
+        cycle_length.append(0)
+        for idx, head in enumerate(head_sequence_result[cycle]):
+            point_index = placement[head]
+            if point_index == -1:
+                continue
+            pos = [mount_point_pos[point_index][0] - head * head_interval, mount_point_pos[point_index][1]]
+            angle = mount_point_angle[point_index]
+            cycle_pos_list.append(pos)
+            if prev_pos is not None:
+                if head_sequence_result[cycle][idx - 1] // 2 == head_sequence_result[cycle][idx] // 2:  # 同轴
+                    rotary_angle = prev_angle - angle
+                else:
+                    rotary_angle = 0
+
+                cycle_length[-1] += max(axis_moving_time(prev_pos[0] - pos[0], 0),
+                                        axis_moving_time(prev_pos[1] - pos[1], 1), head_rotary_time(rotary_angle))
+            prev_pos, prev_angle = pos, angle
+
+        cycle_average_pos.append([sum(map(lambda pos: pos[0], cycle_pos_list)) / len(cycle_pos_list),
+                                  sum(map(lambda pos: pos[1], cycle_pos_list)) / len(cycle_pos_list)])
+
+    best_placement_result, best_head_sequence_result = copy.deepcopy(placement_result), copy.deepcopy(
+        head_sequence_result)
+
+    best_cycle_length, best_cycle_average_pos = copy.deepcopy(cycle_length), copy.deepcopy(cycle_average_pos)
+
+    n_runningtime, n_iteration = 10, 0
+    start_time = time.time()
+    with tqdm(total=n_runningtime, leave=False) as pbar:
+        pbar.set_description('swap heuristic process')
+        prev_time = start_time
+        while True:
+            n_iteration += 1
+
+            placement_result, head_sequence_result = copy.deepcopy(best_placement_result), copy.deepcopy(
+                best_head_sequence_result)
+            cycle_length = best_cycle_length.copy()
+            cycle_average_pos = copy.deepcopy(best_cycle_average_pos)
+
+            cycle_index = roulette_wheel_selection(cycle_length)  # 根据周期加权移动距离随机选择周期
+
+            point_dist = []     # 周期内各贴装点距离中心位置的切氏距离
+            for head in head_sequence_result[cycle_index]:
+                point_index = placement_result[cycle_index][head]
+                _delta_x = abs(mount_point_pos[point_index][0] - head * head_interval - cycle_average_pos[cycle_index][0])
+                _delta_y = abs(mount_point_pos[point_index][1] - cycle_average_pos[cycle_index][1])
+                point_dist.append(max(_delta_x, _delta_y))
+
+            # 随机选择一个异常点
+            head_index = head_sequence_result[cycle_index][roulette_wheel_selection(point_dist)]
+            point_index = placement_result[cycle_index][head_index]
+
+            # 找距离该异常点最近的周期
+            min_dist = None
+            chg_cycle_index = -1
+            for idx in range(len(cycle_average_pos)):
+                if idx == cycle_index:
+                    continue
+                dist_ = 0
+                component_type_check = False
+                for head in head_sequence_result[idx]:
+                    dist_ += max(abs(mount_point_pos[placement_result[idx][head]][0] - mount_point_pos[point_index][0]),
+                                 abs(mount_point_pos[placement_result[idx][head]][1] - mount_point_pos[point_index][1]))
+                    if mount_point_part[placement_result[idx][head]] == mount_point_part[point_index]:
+                        component_type_check = True
+
+                if (min_dist is None or dist_ < min_dist) and component_type_check:
+                    min_dist = dist_
+                    chg_cycle_index = idx
+
+            assert chg_cycle_index != -1
+
+            chg_head, min_chg_dist = None, None
+            chg_cycle_point = []
+            for head in head_sequence_result[chg_cycle_index]:
+                index = placement_result[chg_cycle_index][head]
+                chg_cycle_point.append([mount_point_pos[index][0] - head * head_interval, mount_point_pos[index][1]])
+
+            for idx, head in enumerate(head_sequence_result[chg_cycle_index]):
+                chg_cycle_point_cpy = copy.deepcopy(chg_cycle_point)
+                index = placement_result[chg_cycle_index][head]
+                if mount_point_part[index] != mount_point_part[point_index]:
+                    continue
+                chg_cycle_point_cpy[idx][0] = (mount_point_pos[index][0]) - head * head_interval
+
+                chg_dist = 0
+                aver_chg_pos = [sum(map(lambda x: x[0], chg_cycle_point_cpy)) / len(chg_cycle_point_cpy),
+                                sum(map(lambda x: x[1], chg_cycle_point_cpy)) / len(chg_cycle_point_cpy)]
+
+                for pos in chg_cycle_point_cpy:
+                    chg_dist += max(abs(aver_chg_pos[0] - pos[0]), abs(aver_chg_pos[1] - pos[1]))
+
+                # 更换后各点距离中心更近
+                if min_chg_dist is None or chg_dist < min_chg_dist:
+                    chg_head = head
+                    min_chg_dist = chg_dist
+
+            assert chg_head is not None
+
+            # === 第一轮，变更周期chg_cycle_index的贴装点重排 ===
+            chg_placement_res = placement_result[chg_cycle_index].copy()
+            chg_placement_res[chg_head] = point_index
+
+            cycle_point_list = defaultdict(list)
+            for head, point in enumerate(chg_placement_res):
+                if point == -1:
+                    continue
+                cycle_point_list[mount_point_part[point]].append(point)
+
+            for key, point_list in cycle_point_list.items():
+                cycle_point_list[key] = sorted(point_list, key=lambda p: mount_point_pos[p][0])
+
+            chg_placement_res, chg_point_assign_res = [], [[0, 0]] * max_head_index
+            chg_angle_res = [0] * max_head_index
+            for head, point_index in enumerate(placement_result[chg_cycle_index]):
+                if point_index == -1:
+                    chg_placement_res.append(-1)
+                else:
+                    part = mount_point_part[point_index]
+                    chg_placement_res.append(cycle_point_list[part][0])
+                    chg_point_assign_res[head] = mount_point_pos[cycle_point_list[part][0]].copy()
+                    chg_angle_res[head] = mount_point_angle[cycle_point_list[part][0]]
+                    cycle_point_list[part].pop(0)
+
+            chg_place_moving, chg_head_res = dynamic_programming_cycle_path(chg_placement_res, chg_point_assign_res, chg_angle_res)
+
+            # === 第二轮，原始周期cycle_index的贴装点重排 ===
+            placement_res = placement_result[cycle_index].copy()
+            placement_res[head_index] = placement_result[chg_cycle_index][chg_head]
+
+            for point in placement_res:
+                if point == -1:
+                    continue
+                cycle_point_list[mount_point_part[point]].append(point)
+
+            for key, point_list in cycle_point_list.items():
+                cycle_point_list[key] = sorted(point_list, key=lambda p: mount_point_pos[p][0])
+
+            placement_res, point_assign_res = [], [[0, 0]] * max_head_index
+            angle_assign_res = [0] * max_head_index
+            for head, point_index in enumerate(placement_result[cycle_index]):
+                if point_index == -1:
+                    placement_res.append(-1)
+                else:
+                    part = mount_point_part[point_index]
+                    placement_res.append(cycle_point_list[part][0])
+                    point_assign_res[head] = mount_point_pos[cycle_point_list[part][0]].copy()
+                    angle_assign_res[head] = mount_point_angle[cycle_point_list[part][0]]
+                    cycle_point_list[part].pop(0)
+
+            place_moving, place_head_res = dynamic_programming_cycle_path(placement_res, point_assign_res, angle_assign_res)
+
+            # 更新贴装顺序分配结果
+            placement_result[cycle_index], head_sequence_result[cycle_index] = placement_res, place_head_res
+            placement_result[chg_cycle_index], head_sequence_result[chg_cycle_index] = chg_placement_res, chg_head_res
+
+            # 更新移动路径
+            cycle_length[cycle_index], cycle_length[chg_cycle_index] = place_moving, chg_place_moving
+
+            # 更新平均坐标和最大偏离点索引
+            point_list, point_index_list = [], []
+            for head in head_sequence_result[cycle_index]:
+                point_index_list.append(placement_result[cycle_index][head])
+                point_pos = mount_point_pos[point_index_list[-1]].copy()
+                point_pos[0] -= head * head_interval
+                point_list.append(point_pos)
+
+            cycle_average_pos[cycle_index] = [sum(map(lambda x: x[0], point_list)) / len(point_list),
+                                              sum(map(lambda x: x[1], point_list)) / len(point_list)]
+
+            point_list, point_index_list = [], []
+            for head in head_sequence_result[chg_cycle_index]:
+                point_index_list.append(placement_result[chg_cycle_index][head])
+                point_pos = mount_point_pos[point_index_list[-1]].copy()
+                point_pos[0] -= head * head_interval
+                point_list.append(point_pos)
+
+            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)]
+
+            if sum(cycle_length) < sum(best_cycle_length):
+                best_cycle_length = cycle_length.copy()
+                best_cycle_average_pos = copy.deepcopy(cycle_average_pos)
+                best_placement_result, best_head_sequence_result = copy.deepcopy(placement_result), copy.deepcopy(
+                    head_sequence_result)
+
+            cur_time = time.time()
+            if cur_time - start_time > n_runningtime:
+                break
+
+            pbar.update(cur_time - prev_time)
+            prev_time = cur_time
+
+    # print("number of iteration: ", n_iteration)
+    return best_placement_result, best_head_sequence_result

base_optimizer/result_analysis.py

@@ -1,6 +1,59 @@
 from base_optimizer.optimizer_common import *
 
 
+def convert_pcbdata_to_result(pcb_data, component_data):
+    component_result, cycle_result, feeder_slot_result = [], [], []
+    placement_result, head_sequence_result = [], []
+
+    assigned_part = [-1 for _ in range(max_head_index)]
+    assigned_slot = [-1 for _ in range(max_head_index)]
+    assigned_point = [-1 for _ in range(max_head_index)]
+    assigned_sequence = []
+
+    point_num = len(pcb_data)           # total mount points num
+    for point_cnt in range(point_num + 1):
+
+        cycle_start = 1 if point_cnt == point_num else pcb_data.loc[point_cnt, 'cs']
+        if (cycle_start and point_cnt != 0) or  -1 not in assigned_part:
+
+            if len(component_result) != 0 and component_result[-1] == assigned_part:
+                cycle_result[-1] += 1
+            else:
+                component_result.append(assigned_part)
+                feeder_slot_result.append(assigned_slot)
+                cycle_result.append(1)
+
+            # assigned_sequence = list(reversed(assigned_sequence))       # Samsung拾取顺序相反
+
+            placement_result.append(assigned_point)
+            head_sequence_result.append(assigned_sequence)
+
+            assigned_part = [-1 for _ in range(max_head_index)]
+            assigned_slot = [-1 for _ in range(max_head_index)]
+            assigned_point = [-1 for _ in range(max_head_index)]
+            assigned_sequence = []
+
+            if point_cnt == point_num:
+                break
+
+        slot = pcb_data.loc[point_cnt, 'fdr'].split(' ')[0]
+        if slot == 'A':
+            slot, part = 0, pcb_data.loc[point_cnt].part
+        else:
+
+            slot, part = int(slot[1:]), pcb_data.loc[point_cnt].fdr.split(' ', 1)[1]
+        head = pcb_data.loc[point_cnt].hd - 1
+
+        part_index = component_data[component_data.part == part].index.tolist()[0]
+
+        assigned_part[head] = part_index
+        assigned_slot[head] = slot
+        assigned_point[head] = point_cnt
+        assigned_sequence.append(head)
+
+    return component_result, cycle_result, feeder_slot_result, placement_result, head_sequence_result
+
+
 # 绘制各周期从供料器周期拾取的元件位置
 def pickup_cycle_schematic(feeder_slot_result, cycle_result):
     plt.rcParams['font.sans-serif'] = ['KaiTi']  # 指定默认字体
@@ -329,32 +382,6 @@ def output_optimize_result(file_name, method, component_data, pcb_data, feeder_d
     output_data.to_excel('result/' + file_name, sheet_name='tb1', float_format='%.3f', na_rep='')
 
 
-def component_assign_evaluate(component_data, component_result, cycle_result, feeder_slot_result) -> float:
-    nozzle_change_counter = 0
-    for head in range(max_head_index):
-        nozzle = ''
-        for cycle in range(len(component_result)):
-            component_index = component_result[cycle][head]
-            if component_index == -1:
-                continue
-
-            if cycle != 0 and nozzle != component_data.loc[component_index, 'nz']:
-                nozzle_change_counter += 1
-            nozzle = component_data.loc[component_index, 'nz']
-
-    gang_pick_counter = 0
-    for cycle, feeder_slot in enumerate(feeder_slot_result):
-        pick_slot = defaultdict(int)
-        for head, slot in enumerate(feeder_slot):
-            if slot == -1:
-                continue
-            pick_slot[slot - head * interval_ratio] += 1
-        for _ in pick_slot.values():
-            gang_pick_counter += cycle_result[cycle]
-
-    return sum(cycle_result) + e_nz_change * nozzle_change_counter + e_gang_pick * gang_pick_counter
-
-
 def optimization_assign_result(component_data, pcb_data, component_result, cycle_result, feeder_slot_result,
                                nozzle_hinter=False, component_hinter=False, feeder_hinter=False):
     if nozzle_hinter:
@@ -549,9 +576,9 @@ def placement_time_estimate(component_data, pcb_data, component_result, cycle_re
                 index = placement_result[cycle][head]
                 if index == -1:
                     continue
-                mount_pos.append([pcb_data.loc[index]['x'] - head * head_interval + stopper_pos[0],
-                                  pcb_data.loc[index]['y'] + stopper_pos[1]])
-                mount_angle.append(pcb_data.loc[index]['r'])
+                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'])
 
             # 单独计算贴装路径
             for cntPoints in range(len(mount_pos) - 1):

dataloader.py

@@ -1,5 +1,3 @@
-import random
-
 from base_optimizer.optimizer_common import *
 
 
@@ -51,9 +49,9 @@ def load_data(filename: str, default_feeder_limit=1, load_cp_data=True, load_fee
                 # warnings.warn(warning_info, UserWarning)
         part_index = component_data[component_data['part'] == part].index.tolist()[0]
         part_feeder_assign[part].add(slot)
-        component_data.loc[part_index]['points'] += 1
+        component_data.loc[part_index, 'points'] += 1
 
-        if nozzle != 'A' and component_data.loc[part_index]['nz'] != nozzle:
+        if nozzle != 'A' and component_data.loc[part_index, 'nz'] != nozzle:
             warning_info = 'the nozzle type of component ' + part + ' is not consistent with the pcb data'
             warnings.warn(warning_info, UserWarning)
 

generator.py

@@ -0,0 +1,156 @@
+import pandas as pd
+
+from base_optimizer.optimizer_common import *
+
+
+class DataMgr:
+    def __init__(self):
+        self.min_placement_points = 100
+        self.max_placement_points = 800
+
+        self.max_component_types = 50
+        self.default_feeder_limit = 1
+        self.nozzle_type_list = ['CN065', 'CN140', 'CN220', 'CN040']
+
+        self.x_range = [50, 100, 150, 200, 300, 400, 500]
+        self.y_range = [50, 100, 150, 200, 300, 400, 500]
+
+        self.counter = 0
+        self.pre_file = None
+
+    def generator(self):
+        boundary = [random.choice(self.x_range), random.choice(self.y_range)]
+        total_points = random.randint(self.min_placement_points, self.max_placement_points)
+
+        # determine the nozzle type of component
+        component_list = defaultdict(str)
+        for cp_idx in range(min(random.randint(1, self.max_component_types), total_points)):
+            component_list['C' + str(cp_idx)] = random.choice(self.nozzle_type_list)
+
+        step_col = ["ref", "x", "y", "z", "r", "part", "desc", "fdr", "nz", "hd", "cs", "cy", "sk", "bl", "ar", "pl", "lv"]
+        pcb_data = pd.DataFrame(columns=step_col)
+
+        for idx in range(total_points):
+            part = random.choice(list(component_list.keys()))
+            nozzle = component_list[part]
+
+            pos_x, pos_y = np.random.uniform(0, boundary[0]), np.random.uniform(0, boundary[1])
+            pcb_data = pd.concat([pcb_data, pd.DataFrame([['R' + str(idx), -pos_x, pos_y,
+                                                       0.000, 0.000, part, '', 'A', '1-0 ' + nozzle, 1, 1, 1, 0,
+                                                       1, 1, 1, 'L0']], columns=pcb_data.columns)], ignore_index=True)
+
+        part_col = ["part", "desc", "fdr", "nz", 'camera', 'group', 'feeder-limit', 'points']
+        component_data = pd.DataFrame(columns=part_col)
+
+        for _, data in pcb_data.iterrows():
+            part, nozzle = data.part, data.nz.split(' ')[1]
+            if part not in component_data['part'].values:
+                component_data = pd.concat([component_data, pd.DataFrame(
+                    [part, '', 'SM8', nozzle, '飞行相机1', 'CHIP-Rect', self.default_feeder_limit, 0], index=part_col).T],
+                                           ignore_index=True)
+
+            part_index = component_data[component_data['part'] == part].index.tolist()[0]
+            component_data.loc[part_index, 'points'] += 1
+        self.counter += 1
+        return pcb_data, component_data
+
+    def recorder(self, file_path, info: OptInfo, pcb_data, component_data):
+        lineinfo = '{:.6f}'.format(info.placement_time) + '\t' + str(info.cycle_counter) + '\t' + str(
+            info.nozzle_change_counter) + '\t' + str(info.pickup_counter) + '\t' + '{:.3f}'.format(
+            info.pickup_movement) + '\t' + '{:.3f}'.format(info.placement_movement)
+
+        lineinfo += '\t' + '{:.3f}'.format(pcb_data['x'].max() - pcb_data['x'].min()) + '\t' + '{:.3f}'.format(
+            pcb_data['y'].max() - pcb_data['y'].min())
+
+        point_counter, component_counter = 0, 0
+        nozzle_type = set()
+        for _, data in component_data.iterrows():
+            if data.points == 0:
+                continue
+            nozzle_type.add(data.nz)
+            point_counter += data.points
+            component_counter += 1
+
+        lineinfo += '\t' + str(point_counter) + '\t' + str(component_counter) + '\t' + str(len(nozzle_type))
+
+        for _, data in component_data.iterrows():
+            lineinfo += '\t' + data.part + '\t' + data.nz + '\t' + str(data.points)
+        lineinfo += '\n'
+
+        with open(file_path, 'a') as f:
+            f.write(lineinfo)
+        f.close()
+
+    def saver(self, file_path: str, pcb_data):
+        lineinfo = ''
+        for _, data in pcb_data.iterrows():
+            lineinfo += '\t' + '{:.3f}'.format(data.x) + '\t' + '{:.3f}'.format(
+                data.y) + '\t0.000\t0.000\t' + data.part + '\t\tA\t' + data.nz + '\t1\t1\t1\t1\t1\t1\t1\tN\tL0\n'
+        pos = file_path.find('.')
+        file_path = file_path[:pos] + '-' + str(self.counter) + file_path[pos:]
+        with open(file_path, 'w') as f:
+            f.write(lineinfo)
+        f.close()
+        self.pre_file = file_path
+
+    def remover(self):
+        if self.pre_file is not None:
+            os.remove(self.pre_file)
+            self.pre_file = None
+
+    def encode(self, cp_points: defaultdict[str], cp_nozzle: defaultdict[int], width, height):
+        cp2nz = defaultdict(int)
+        for idx, nozzle in enumerate(self.nozzle_type_list):
+            cp2nz[nozzle] = idx
+
+        total_points = sum(points for points in cp_points.values())
+        total_component_types, total_nozzle_types = len(cp_points.keys()), len(set(cp_nozzle.values()))
+        data = [total_points, total_component_types, total_nozzle_types]
+        data.extend([width, height])
+
+        for component, points in cp_points.items():
+            nozzle = cp_nozzle[component]
+
+            data_slice = [0 for _ in range(len(self.nozzle_type_list))]
+            data_slice[cp2nz[nozzle]] = points
+            data.extend(data_slice)
+
+        for _ in range(self.max_component_types - total_component_types):
+            data.extend([0 for _ in range(len(self.nozzle_type_list))])
+
+        return data
+
+    def loader(self, file_path):
+        train_data, time_data = [], []
+        cycle_data, nozzle_change_data, pickup_data, movement_data, point_data = [], [], [], [], []
+        pcb_width, pcb_height = [], []
+        with open(file_path, 'r') as file:
+            line = file.readline()
+            while line:
+                items = line.split('\t')
+                total_points, total_component_types = int(items[8]), int(items[9])
+
+                cycle_data.append(float(items[1]))
+                nozzle_change_data.append(float(items[2]))
+                pickup_data.append(float(items[3]))
+                movement_data.append(float(items[4]))
+                point_data.append(total_points)
+
+                # assembly time data
+                time_data.append(float(items[0]))
+
+                cp_points, cp_nozzle = defaultdict(int), defaultdict(str)
+                for cp_counter in range(total_component_types):
+                    component_type, nozzle_type = items[11 + cp_counter * 3], items[12 + cp_counter * 3]
+                    points = int(items[13 + cp_counter * 3])
+
+                    cp_points[component_type], cp_nozzle[component_type] = points, nozzle_type
+
+                train_data.append(self.encode(cp_points, cp_nozzle, float(items[6]), float(items[7])))
+                line = file.readline()
+
+        return train_data, time_data, cycle_data, nozzle_change_data, pickup_data, movement_data, point_data
+
+    def get_feature(self):
+        return self.max_component_types * len(self.nozzle_type_list) + 5
+

optimizer.py

@@ -1,18 +1,8 @@
-import copy
-import math
-
-import matplotlib.pyplot as plt
-import pandas as pd
-
-from base_optimizer.optimizer_scanbased import *
-from base_optimizer.optimizer_celldivision import *
-from base_optimizer.optimizer_hybridgenetic import *
-from base_optimizer.optimizer_feederpriority import *
 from dataloader import *
-
 from optimizer_genetic import *
 from optimizer_heuristic import *
 from optimizer_reconfiguration import *
+from base_optimizer.optimizer_interface import *
 
 
 def deviation(data):
@@ -23,16 +13,20 @@ def deviation(data):
     return variance / len(data)
 
 
-def optimizer(pcb_data, component_data, assembly_line_optimizer, single_machine_optimizer):
-    # todo: 由于吸嘴更换更因素的存在，在处理PCB8数据时，遗传算法因在负载均衡过程中对这一因素进行了考虑，性能更优
-    # assignment_result = assemblyline_optimizer_heuristic(pcb_data, component_data)
-    # assignment_result = assemblyline_optimizer_genetic(pcb_data, component_data)
-    assignment_result = reconfiguration_optimizer(pcb_data, component_data)
+def optimizer(pcb_data, component_data, line_optimizer, machine_optimizer, machine_number):
+    if line_optimizer == "heuristic":
+        assignment_result = assemblyline_optimizer_heuristic(pcb_data, component_data, machine_number)
+    elif line_optimizer == "genetic":
+        assignment_result = assemblyline_optimizer_genetic(pcb_data, component_data, machine_number)
+    elif line_optimizer == "reconfiguration":
+        assignment_result = reconfiguration_optimizer(pcb_data, component_data, machine_number)
+    else:
+        return
 
     assignment_result_cpy = copy.deepcopy(assignment_result)
     placement_points, placement_time = [], []
     partial_pcb_data, partial_component_data = defaultdict(pd.DataFrame), defaultdict(pd.DataFrame)
-    for machine_index in range(max_machine_index):
+    for machine_index in range(machine_number):
         partial_pcb_data[machine_index] = pd.DataFrame(columns=pcb_data.columns)
         partial_component_data[machine_index] = component_data.copy(deep=True)
         placement_points.append(sum(assignment_result[machine_index]))
@@ -42,29 +36,29 @@ def optimizer(pcb_data, component_data, assembly_line_optimizer, single_machine_
     # === averagely assign available feeder ===
     for part_index, data in component_data.iterrows():
         feeder_limit = data['feeder-limit']
-        feeder_points = [assignment_result[machine_index][part_index] for machine_index in range(max_machine_index)]
+        feeder_points = [assignment_result[machine_index][part_index] for machine_index in range(machine_number)]
 
-        for machine_index in range(max_machine_index):
+        for machine_index in range(machine_number):
             if feeder_points[machine_index] == 0:
                 continue
 
             arg_feeder = max(math.floor(feeder_points[machine_index] / sum(feeder_points) * data['feeder-limit']), 1)
 
-            partial_component_data[machine_index].loc[part_index]['feeder-limit'] = arg_feeder
+            partial_component_data[machine_index].loc[part_index, 'feeder-limit'] = arg_feeder
             feeder_limit -= arg_feeder
 
-        for machine_index in range(max_machine_index):
+        for machine_index in range(machine_number):
             if feeder_limit <= 0:
                 break
 
             if feeder_points[machine_index] == 0:
                 continue
-            partial_component_data[machine_index].loc[part_index]['feeder-limit'] += 1
+            partial_component_data[machine_index].loc[part_index, 'feeder-limit'] += 1
             feeder_limit -= 1
 
-        for machine_index in range(max_machine_index):
+        for machine_index in range(machine_number):
             if feeder_points[machine_index] > 0:
-                assert partial_component_data[machine_index].loc[part_index]['feeder-limit'] > 0
+                assert partial_component_data[machine_index].loc[part_index, 'feeder-limit'] > 0
 
     # === assign placements ===
     component_machine_index = [0 for _ in range(len(component_data))]
@@ -102,7 +96,7 @@ def optimizer(pcb_data, component_data, assembly_line_optimizer, single_machine_
                     assign_part_index.append(idx_)
 
             variance = deviation(assign_part_point)
-            while start_index != end_index:
+            while start_index <= end_index:
                 part_info_index = assign_part_index[np.argmax(assign_part_point)]
 
                 if part_info[part_info_index][2] < part_info[part_info_index][3]:   # 供料器数目上限的限制
@@ -122,7 +116,7 @@ def optimizer(pcb_data, component_data, assembly_line_optimizer, single_machine_
                         break
 
                     variance = new_variance
-                    assign_part_index, assign_part_point = new_assign_part_index, new_assign_part_point
+                    assign_part_index, assign_part_point = new_assign_part_index.copy(), new_assign_part_point.copy()
                 else:
                     break
 
@@ -132,108 +126,39 @@ def optimizer(pcb_data, component_data, assembly_line_optimizer, single_machine_
         # update available feeder number
         max_avl_feeder = max(part_info, key=lambda x: x[2])[2]
         for info in part_info:
-            partial_component_data[machine_index].loc[info[0]]['feeder-limit'] = math.ceil(info[2] / max_avl_feeder)
+            partial_component_data[machine_index].loc[info[0], 'feeder-limit'] = math.ceil(info[2] / max_avl_feeder)
 
         placement_time.append(base_optimizer(machine_index + 1, data, partial_component_data[machine_index],
                                              feeder_data=pd.DataFrame(columns=['slot', 'part', 'arg']),
-                                             method=single_machine_optimizer, hinter=True))
+                                             method=machine_optimizer, hinter=True).placement_time)
 
-    average_time, standard_deviation_time = sum(placement_time) / max_machine_index, 0
-    for machine_index in range(max_machine_index):
+    average_time, standard_deviation_time = sum(placement_time) / machine_number, 0
+    for machine_index in range(machine_number):
+        total_component_types = 0
+        for points in assignment_result_cpy[machine_index]:
+            if points:
+                total_component_types += 1
         print('assembly time for machine ' + str(machine_index + 1) + ': ' + str(
-            placement_time[machine_index]) + ' s, ' + 'total placements: ' + str(placement_points[machine_index]))
+            placement_time[machine_index]) + ' s, ' + 'total placements: ' + str(placement_points[machine_index])
+              + ', total component types: ' + str(total_component_types))
         standard_deviation_time += pow(placement_time[machine_index] - average_time, 2)
-    standard_deviation_time /= max_machine_index
+    standard_deviation_time /= machine_number
     standard_deviation_time = math.sqrt(standard_deviation_time)
 
     print('finial assembly time: ' + str(max(placement_time)) + 's, standard deviation: ' + str(standard_deviation_time))
 
 
-# todo: 不同类型元件的组装时间差异
-def base_optimizer(machine_index, pcb_data, component_data, feeder_data=None, method='', hinter=False):
-
-    if method == 'cell_division':  # 基于元胞分裂的遗传算法
-        component_result, cycle_result, feeder_slot_result = optimizer_celldivision(pcb_data, component_data,
-                                                                                    hinter=False)
-        placement_result, head_sequence = greedy_placement_route_generation(component_data, pcb_data, component_result,
-                                                                            cycle_result, feeder_slot_result)
-    elif method == 'feeder_scan':  # 基于基座扫描的供料器优先算法
-        # 第1步：分配供料器位置
-        nozzle_pattern = feeder_allocate(component_data, pcb_data, feeder_data, figure=False)
-        # 第2步：扫描供料器基座，确定元件拾取的先后顺序
-        component_result, cycle_result, feeder_slot_result = feeder_base_scan(component_data, pcb_data, feeder_data,
-                                                                              nozzle_pattern)
-
-        # 第3步：贴装路径规划
-        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_for_route_generation(component_data, pcb_data, component_result,
-        #                                                                    cycle_result, feeder_slot_result)
-
-    elif method == 'hybrid_genetic':  # 基于拾取组的混合遗传算法
-        component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_hybrid_genetic(
-            pcb_data, component_data, hinter=False)
-
-    elif method == 'aggregation':  # 基于batch-level的整数规划 + 启发式算法
-        component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_aggregation(
-            component_data, pcb_data)
-    elif method == 'genetic_scanning':
-        component_result, cycle_result, feeder_slot_result, placement_result, head_sequence = optimizer_genetic_scanning(
-            component_data, pcb_data, hinter=False)
-    else:
-        raise 'method is not existed'
-
-    if hinter:
-        optimization_assign_result(component_data, pcb_data, component_result, cycle_result, feeder_slot_result,
-                                   nozzle_hinter=True, component_hinter=False, feeder_hinter=False)
-
-        print('----- Placement machine ' + str(machine_index) + ' ----- ')
-        print('-Cycle counter: {}'.format(sum(cycle_result)))
-
-        total_nozzle_change_counter, total_pick_counter = 0, 0
-        total_pick_movement = 0
-        assigned_nozzle = ['' if idx == -1 else component_data.loc[idx]['nz'] for idx in component_result[0]]
-
-        for cycle in range(len(cycle_result)):
-            pick_slot = set()
-            for head in range(max_head_index):
-                if (idx := component_result[cycle][head]) == -1:
-                    continue
-
-                nozzle = component_data.loc[idx]['nz']
-                if nozzle != assigned_nozzle[head]:
-                    if assigned_nozzle[head] != '':
-                        total_nozzle_change_counter += 1
-                    assigned_nozzle[head] = nozzle
-
-                pick_slot.add(feeder_slot_result[cycle][head] - head * interval_ratio)
-            total_pick_counter += len(pick_slot) * cycle_result[cycle]
-
-            pick_slot = list(pick_slot)
-            pick_slot.sort()
-            for idx in range(len(pick_slot) - 1):
-                total_pick_movement += abs(pick_slot[idx+1] - pick_slot[idx])
-
-        print('-Nozzle change counter: {}'.format(total_nozzle_change_counter))
-        print('-Pick operation counter: {}'.format(total_pick_counter))
-        print('-Pick movement: {}'.format(total_pick_movement))
-        print('------------------------------ ')
-
-    # 估算贴装用时
-    return placement_time_estimate(component_data, pcb_data, component_result, cycle_result, feeder_slot_result,
-                                   placement_result, head_sequence, hinter=False)
-
-
 @timer_wrapper
 def main():
-    # warnings.simplefilter('ignore')
+    warnings.simplefilter(action='ignore', category=FutureWarning)
     # 参数解析
     parser = argparse.ArgumentParser(description='assembly line optimizer implementation')
     parser.add_argument('--filename', default='PCB.txt', type=str, help='load pcb data')
     parser.add_argument('--auto_register', default=1, type=int, help='register the component according the pcb data')
-    parser.add_argument('--base_optimizer', default='feeder_scan', type=str, help='base optimizer for single machine')
-    parser.add_argument('--assembly_optimizer', default='heuristic', type=str, help='optimizer for PCB Assembly Line')
-    parser.add_argument('--feeder_limit', default=2, type=int,
+    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_scan', type=str, help='optimizer for single machine')
+    parser.add_argument('--line_optimizer', default='genetic', 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')
     params = parser.parse_args()
 
@@ -243,9 +168,8 @@ def main():
 
     # 加载PCB数据
     pcb_data, component_data, _ = load_data(params.filename, default_feeder_limit=params.feeder_limit,
-                                            cp_auto_register=params.auto_register)  # 加载PCB数据
-
-    optimizer(pcb_data, component_data, params.assembly_optimizer, params.base_optimizer)
+                                            cp_auto_register=params.auto_register, load_feeder_data=False)  # 加载PCB数据
+    optimizer(pcb_data, component_data, params.line_optimizer, params.machine_optimizer, params.machine_number)
 
 
 if __name__ == '__main__':

optimizer_genetic.py

@@ -1,12 +1,9 @@
 # implementation of <<An integrated allocation method for the PCB assembly line balancing problem with nozzle changes>>
-import copy
-
-import matplotlib.pyplot as plt
-
 from base_optimizer.optimizer_common import *
+from optimizer_hyperheuristic import *
 
 
-def selective_initialization(component_points, component_feeders, population_size):
+def selective_initialization(component_points, component_feeders, population_size, machine_number):
     population = []    # population initialization
     for _ in range(population_size):
         individual = []
@@ -14,14 +11,14 @@ def selective_initialization(component_points, component_feeders, population_siz
             if points == 0:
                 continue
             # 可用机器数
-            avl_machine_num = random.randint(1, min(max_machine_index, component_feeders[part_index], points))
+            avl_machine_num = random.randint(1, min(machine_number, component_feeders[part_index], points))
 
             selective_possibility = []
             for p in range(1, avl_machine_num + 1):
                 selective_possibility.append(pow(2, avl_machine_num - p + 1))
 
             sel_machine_num = random_selective([p + 1 for p in range(avl_machine_num)], selective_possibility)  # 选择的机器数
-            sel_machine_set = random.sample([p for p in range(max_machine_index)], sel_machine_num)
+            sel_machine_set = random.sample([p for p in range(machine_number)], sel_machine_num)
 
             sel_machine_points = [1 for _ in range(sel_machine_num)]
             for p in range(sel_machine_num - 1):
@@ -34,7 +31,7 @@ def selective_initialization(component_points, component_feeders, population_siz
                 sel_machine_points[-1] += (points - sum(sel_machine_points))
 
             # code component allocation into chromosome
-            for p in range(max_machine_index):
+            for p in range(machine_number):
                 if p in sel_machine_set:
                     individual += [0 for _ in range(sel_machine_points[0])]
                     sel_machine_points.pop(0)
@@ -45,7 +42,7 @@ def selective_initialization(component_points, component_feeders, population_siz
     return population
 
 
-def selective_crossover(component_points, component_feeders, mother, father, non_decelerating=True):
+def selective_crossover(component_points, component_feeders, mother, father, machine_number, non_decelerating=True):
     assert len(mother) == len(father)
 
     offspring1, offspring2 = mother.copy(), father.copy()
@@ -56,24 +53,24 @@ def selective_crossover(component_points, component_feeders, mother, father, non
         one_counter = 0
 
         idx_, mother_cut_line, father_cut_line = 0, [-1], [-1]
-        for idx_, gene in enumerate(mother[idx: idx + points + max_machine_index - 1]):
+        for idx_, gene in enumerate(mother[idx: idx + points + machine_number - 1]):
             if gene:
                 mother_cut_line.append(idx_)
         mother_cut_line.append(idx_ + 1)
 
-        for idx_, gene in enumerate(father[idx: idx + points + max_machine_index - 1]):
+        for idx_, gene in enumerate(father[idx: idx + points + machine_number - 1]):
             if gene:
                 father_cut_line.append(idx_)
         father_cut_line.append(idx_ + 1)
 
-        for offset in range(points + max_machine_index - 1):
+        for offset in range(points + machine_number - 1):
             if mother[idx + offset] == 1:
                 one_counter += 1
             if father[idx + offset] == 1:
                 one_counter -= 1
 
             # first constraint: the total number of '1's (the number of partitions) in the chromosome is unchanged
-            if one_counter != 0 or offset == 0 or offset == points + max_machine_index - 2:
+            if one_counter != 0 or offset == 0 or offset == points + machine_number - 2:
                 continue
 
             # the selected cut-line should guarantee there are the same or a larger number unassigned machine
@@ -89,13 +86,14 @@ def selective_crossover(component_points, component_feeders, mother, father, non
                 n_new += 1
 
             # second constraint: non_decelerating or accelerating crossover
+            # non_decelerating or accelerating means that the number of machine without workload is increased
             if n_new < n_bro or (n_new == n_bro and not non_decelerating):
                 continue
 
             # third constraint (customized constraint):
             # no more than the maximum number of available machine for each component type
             new_mother_cut_line, new_father_cut_line = [], []
-            for idx_ in range(max_machine_index + 1):
+            for idx_ in range(machine_number + 1):
                 if mother_cut_line[idx_] <= offset:
                     new_mother_cut_line.append(mother_cut_line[idx_])
                 else:
@@ -110,11 +108,11 @@ def selective_crossover(component_points, component_feeders, mother, father, non
             sorted(new_father_cut_line, reverse=False)
             n_mother_machine, n_father_machine = 0, 0
 
-            for idx_ in range(max_machine_index):
-                if new_mother_cut_line[idx_ + 1] - new_mother_cut_line[idx_]:
+            for idx_ in range(machine_number):
+                if new_mother_cut_line[idx_ + 1] - new_mother_cut_line[idx_] > 1:
                     n_mother_machine += 1
 
-                if new_father_cut_line[idx_ + 1] - new_father_cut_line[idx_]:
+                if new_father_cut_line[idx_ + 1] - new_father_cut_line[idx_] > 1:
                     n_father_machine += 1
 
             if n_mother_machine > component_feeders[part_index] or n_father_machine > component_feeders[part_index]:
@@ -122,7 +120,7 @@ def selective_crossover(component_points, component_feeders, mother, father, non
 
             feasible_cut_line.append(idx + offset)
 
-        idx += (points + max_machine_index - 1)
+        idx += (points + machine_number - 1)
 
     if len(feasible_cut_line) == 0:
         return offspring1, offspring2
@@ -133,14 +131,14 @@ def selective_crossover(component_points, component_feeders, mother, father, non
     return offspring1, offspring2
 
 
-def cal_individual_val(component_points, component_nozzle, individual):
+def cal_individual_val(component_points, component_feeders, component_nozzle, machine_number, individual, data_mgr, net):
     idx, objective_val = 0, []
-    machine_component_points = [[] for _ in range(max_machine_index)]
+    machine_component_points = [[] for _ in range(machine_number)]
     nozzle_component_points = defaultdict(list)
 
     # decode the component allocation
     for comp_idx, points in component_points:
-        component_gene = individual[idx: idx + points + max_machine_index - 1]
+        component_gene = individual[idx: idx + points + machine_number - 1]
         machine_idx, component_counter = 0, 0
         for gene in component_gene:
             if gene:
@@ -150,14 +148,31 @@ def cal_individual_val(component_points, component_nozzle, individual):
             else:
                 component_counter += 1
         machine_component_points[-1].append(component_counter)
-        idx += (points + max_machine_index - 1)
+        idx += (points + machine_number - 1)
 
         nozzle_component_points[component_nozzle[comp_idx]] = [0] * len(component_points)   # 初始化元件-吸嘴点数列表
 
+    # ======== 新加的开始 ========
+    for machine_idx in range(machine_number):
+        cp_points, cp_nozzle = defaultdict(int), defaultdict(str)
+        for comp_idx, _ in component_points:
+            if machine_component_points[machine_idx][comp_idx] == 0:
+                continue
+            cp_points['C' + str(comp_idx)] = machine_component_points[machine_idx][comp_idx]
+            cp_nozzle['C' + str(comp_idx)] = component_nozzle[comp_idx]
+
+        encoding = np.array(data_mgr.encode(cp_points, cp_nozzle, 45, 150))
+        encoding = torch.from_numpy(encoding.reshape((-1, np.shape(encoding)[0]))).float().to("cuda")
+        # pred_time = net(encoding)[0, 0].item()
+        # objective_val.append(pred_time * sum(points for points in cp_points.values()))
+        objective_val.append(net(encoding)[0, 0].item())
+
+    return objective_val, machine_component_points
+    # ======== 新加的结束（以下内容弃用） =====
     for comp_idx, points in component_points:
         nozzle_component_points[component_nozzle[comp_idx]][comp_idx] = points
 
-    for machine_idx in range(max_machine_index):
+    for machine_idx in range(machine_number):
         nozzle_points = defaultdict(int)
         for idx, nozzle in component_nozzle.items():
             if component_points[idx] == 0:
@@ -236,20 +251,39 @@ def cal_individual_val(component_points, component_nozzle, individual):
         for idx in range(len(heads_placement) // max_head_index):
             wl += heads_placement[idx][1]
         objective_val.append(T_pp * machine_points + T_tr * wl + T_nc * ul + T_pl * pl)
-<<<<<<< HEAD
 
-=======
->>>>>>> 87ddb057cadf152d7af793aa7b8da439dedbe361
     return objective_val, machine_component_points
 
 
-def assemblyline_optimizer_genetic(pcb_data, component_data):
+def individual_convert(component_points, individual):
+    machine_number = len(individual)
+    machine_component_points = [[] for _ in range(machine_number)]
+    idx = 0
+    # decode the component allocation
+    for comp_idx, points in component_points:
+        component_gene = individual[idx: idx + points + machine_number - 1]
+        machine_idx, component_counter = 0, 0
+        for gene in component_gene:
+            if gene:
+                machine_component_points[machine_idx].append(component_counter)
+                machine_idx += 1
+                component_counter = 0
+            else:
+                component_counter += 1
+        machine_component_points[-1].append(component_counter)
+        idx += (points + machine_number - 1)
+
+    return machine_component_points
+
+
+def assemblyline_optimizer_genetic(pcb_data, component_data, machine_number):
     # basic parameter
     # crossover rate & mutation rate: 80% & 10%
     # population size: 200
     # the number of generation: 500
     crossover_rate, mutation_rate = 0.8, 0.1
     population_size, n_generations = 200, 500
+    # population_size, n_generations = 30, 50
 
     # the number of placement points, the number of available feeders, and nozzle type of component respectively
     component_points, component_feeders, component_nozzle = defaultdict(int), defaultdict(int), defaultdict(str)
@@ -262,9 +296,16 @@ def assemblyline_optimizer_genetic(pcb_data, component_data):
         component_nozzle[part_index] = nozzle
 
     component_points = sorted(component_points.items(), key=lambda x: x[0])     # 决定染色体排列顺序
+    data_mgr = DataMgr()
+    device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
+    net = Net(input_size=data_mgr.get_feature(), output_size=1).to(device)
+
+    net.load_state_dict(torch.load('model_state.pth'))
+    # optimizer = torch.optim.Adam(net.parameters(), lr=0.1)
+    # optimizer.load_state_dict(torch.load('optimizer_state.pth'))
 
     # population initialization
-    population = selective_initialization(component_points, component_feeders, population_size)
+    population = selective_initialization(component_points, component_feeders, population_size, machine_number)
     with tqdm(total=n_generations) as pbar:
         pbar.set_description('genetic algorithm process for PCB assembly line balance')
 
@@ -273,7 +314,8 @@ def assemblyline_optimizer_genetic(pcb_data, component_data):
             # calculate fitness value
             pop_val = []
             for individual in population:
-                val, assigned_points = cal_individual_val(component_points, component_nozzle, individual)
+                val, assigned_points = cal_individual_val(component_points, component_feeders, component_nozzle,
+                                                          machine_number, individual, data_mgr, net)
                 pop_val.append(max(val))
 
             select_index = get_top_k_value(pop_val, population_size - len(new_population), reverse=False)
@@ -282,7 +324,8 @@ def assemblyline_optimizer_genetic(pcb_data, component_data):
 
             population += new_population
             for individual in new_population:
-                val, _ = cal_individual_val(component_points, component_nozzle, individual)
+                val, _ = cal_individual_val(component_points, component_feeders, component_nozzle, machine_number,
+                                            individual, data_mgr, net)
                 pop_val.append(max(val))
 
             # min-max convert
@@ -302,13 +345,13 @@ def assemblyline_optimizer_genetic(pcb_data, component_data):
                             break
 
                     offspring1, offspring2 = selective_crossover(component_points, component_feeders,
-                                                                 population[index1], population[index2])
+                                                                 population[index1], population[index2], machine_number)
 
                     if np.random.random() < mutation_rate:
-                        offspring1 = constraint_swap_mutation(component_points, offspring1)
+                        offspring1 = constraint_swap_mutation(component_points, offspring1, machine_number)
 
                     if np.random.random() < mutation_rate:
-                        offspring2 = constraint_swap_mutation(component_points, offspring2)
+                        offspring2 = constraint_swap_mutation(component_points, offspring2, machine_number)
 
                     new_population.append(offspring1)
                     new_population.append(offspring2)
@@ -316,12 +359,13 @@ def assemblyline_optimizer_genetic(pcb_data, component_data):
             pbar.update(1)
 
     best_individual = population[np.argmax(pop_val)]
-    _, assignment_result = cal_individual_val(component_points, component_nozzle, best_individual)
-
+    val, assignment_result = cal_individual_val(component_points, component_feeders, component_nozzle, machine_number,
+                                                best_individual, data_mgr, net)
+    print('final value: ', val)
     # available feeder check
     for part_index, data in component_data.iterrows():
         feeder_limit = data['feeder-limit']
-        for machine_index in range(max_machine_index):
+        for machine_index in range(machine_number):
             if assignment_result[machine_index][part_index]:
                 feeder_limit -= 1
         assert feeder_limit >= 0

optimizer_heuristic.py

@@ -1,13 +1,17 @@
+import copy
 import math
+import random
+
 import numpy as np
 
 from base_optimizer.optimizer_common import *
+from base_optimizer.optimizer_feederpriority import *
+from base_optimizer.result_analysis import *
 
 
 # TODO: nozzle tool available restriction
 # TODO: consider with the PCB placement topology
 def assembly_time_estimator(assignment_points, component_feeders, component_nozzle):
-
     nozzle_heads, nozzle_points = defaultdict(int), defaultdict(int)
     for idx, points in enumerate(assignment_points):
         if points == 0:
@@ -26,27 +30,181 @@ def assembly_time_estimator(assignment_points, component_feeders, component_nozz
         assert max_cycle_nozzle is not None
         nozzle_heads[max_cycle_nozzle] += 1
 
-    n_cycle = max(map(lambda x: math.ceil(nozzle_points[x[0]] / x[1]), nozzle_heads.items()))
+    head_nozzle_assignment, min_cost = None, None
 
-    # calculate the number of simultaneous pickup
-    head_index, nozzle_cycle = 0, [[] for _ in range(max_head_index)]
+    # generate initial nozzle group
+    nozzle_group = []
+    # averagely assign for the same type of nozzles, and generate nozzle group
+    nozzle_points_cpy = copy.deepcopy(nozzle_points)
     for nozzle, heads in nozzle_heads.items():
-        head_index_cpy, points = head_index, nozzle_points[nozzle]
+        points = nozzle_points_cpy[nozzle] // heads
         for _ in range(heads):
-            nozzle_cycle[head_index].append([nozzle, points // heads])
-            head_index += 1
+            nozzle_group.append([nozzle, points])
+        nozzle_points_cpy[nozzle] -= heads * points
 
-        points %= heads
-        while points:
-            nozzle_cycle[head_index_cpy][1] += 1
-            points -= 1
-            head_index_cpy += 1
+    for idx, [nozzle, _] in enumerate(nozzle_group):
+        if nozzle_points_cpy[nozzle]:
+            nozzle_group[idx][1] += 1
+            nozzle_points_cpy[nozzle] -= 1
 
-    # nozzle_cycle_index = [0 for _ in range(max_head_index)]
-    return n_cycle, n_nz_change, n_gang_pick
+    while True:
+        # assign nozzle group to each head
+        nozzle_group.sort(key=lambda x: -x[1])
+
+        tmp_head_nozzle_assignment = []
+        head_total_points = [0 for _ in range(max_head_index)]
+        for idx, nozzle_item in enumerate(nozzle_group):
+            if idx < max_head_index:
+                tmp_head_nozzle_assignment.append([nozzle_item.copy()])
+                head_total_points[idx] += nozzle_item[1]
+            else:
+                min_head = np.argmin(head_total_points)
+                tmp_head_nozzle_assignment[min_head].append(nozzle_item.copy())
+                head_total_points[min_head] += nozzle_item[1]
+
+        cost = t_cycle * max(head_total_points)
+        for head in range(max_head_index):
+            for cycle in range(len(tmp_head_nozzle_assignment[head])):
+                if cycle + 1 == len(tmp_head_nozzle_assignment[head]):
+                    if tmp_head_nozzle_assignment[head][cycle][0] != tmp_head_nozzle_assignment[head][-1][0]:
+                        cost += t_nozzle_change
+                else:
+                    if tmp_head_nozzle_assignment[head][cycle][0] != tmp_head_nozzle_assignment[head][cycle + 1][0]:
+                        cost += t_nozzle_change
+
+        while True:
+            min_head, max_head = np.argmin(head_total_points), np.argmax(head_total_points)
+            min_head_nozzle, max_head_nozzle = tmp_head_nozzle_assignment[min_head][-1][0], \
+                                               tmp_head_nozzle_assignment[max_head][-1][0]
+            if min_head_nozzle == max_head_nozzle:
+                break
+
+            min_head_list, max_head_list = [min_head], [max_head]
+            minmax_head_points = 0
+            for head in range(max_head_index):
+                if head in min_head_list or head in max_head_list:
+                    minmax_head_points += head_total_points[head]
+                    continue
+
+                # the max/min heads with the sum nozzle type
+                if tmp_head_nozzle_assignment[head][-1][0] == tmp_head_nozzle_assignment[min_head][-1][0]:
+                    min_head_list.append(head)
+                    minmax_head_points += head_total_points[head]
+
+                if tmp_head_nozzle_assignment[head][-1][0] == tmp_head_nozzle_assignment[max_head][-1][0]:
+                    max_head_list.append(head)
+                    minmax_head_points += head_total_points[head]
+
+            # todo: restriction of available nozzle
+            # the reduction of cycles is not offset the cost of nozzle change
+            average_points = minmax_head_points // (len(min_head_list) + len(max_head_list))
+            reminder_points = minmax_head_points % (len(min_head_list) + len(max_head_list))
+            max_cycle = average_points + (1 if reminder_points > 0 else 0)
+            for head in range(max_head_index):
+                if head in min_head_list or head in max_head_list:
+                    continue
+                max_cycle = max(max_cycle, head_total_points[head])
+
+            nozzle_change_counter = 0
+            for head in min_head_list:
+                if tmp_head_nozzle_assignment[head][0] == tmp_head_nozzle_assignment[head][-1]:
+                    nozzle_change_counter += 2
+                else:
+                    nozzle_change_counter += 1
+
+            if t_cycle * (max(head_total_points) - max_cycle) < t_nozzle_change * nozzle_change_counter:
+                break
+
+            cost -= t_cycle * (max(head_total_points) - max_cycle) - t_nozzle_change * nozzle_change_counter
+
+            required_points = 0  # 待均摊的贴装点数较多的吸嘴类型
+            for head in min_head_list:
+                points = average_points - head_total_points[head]
+                tmp_head_nozzle_assignment[head].append([max_head_nozzle, points])
+                head_total_points[head] = average_points
+                required_points += points
+
+            for head in max_head_list:
+                tmp_head_nozzle_assignment[head][-1][1] -= required_points // len(max_head_list)
+                head_total_points[head] -= required_points // len(max_head_list)
+
+            required_points -= (required_points // len(max_head_list)) * len(max_head_list)
+
+            for head in max_head_list:
+                if required_points <= 0:
+                    break
+                tmp_head_nozzle_assignment[head][-1][1] -= 1
+                head_total_points[head] -= 1
+                required_points -= 1
+
+        if min_cost is None or cost < min_cost:
+            min_cost = cost
+            head_nozzle_assignment = copy.deepcopy(tmp_head_nozzle_assignment)
+        else:
+            break
+
+        # 在吸嘴组中增加一个吸嘴
+        idx, nozzle = 0, nozzle_group[0][0]
+        for idx, [nozzle_, _] in enumerate(nozzle_group):
+            if nozzle_ != nozzle:
+                break
+
+        average_points, remainder_points = nozzle_points[nozzle] // (idx + 1), nozzle_points[nozzle] % (idx + 1)
+        nozzle_group.append([nozzle, 0])
+        for idx, [nozzle_, _] in enumerate(nozzle_group):
+            if nozzle_ == nozzle:
+                nozzle_group[idx][1] = average_points + (1 if remainder_points > 0 else 0)
+                remainder_points -= 1
+
+    cycle_counter, nozzle_change_counter = 0, 0
+    for head in range(max_head_index):
+        head_cycle_counter = 0
+        for cycle in range(len(head_nozzle_assignment[head])):
+            if cycle + 1 == len(head_nozzle_assignment[head]):
+                if head_nozzle_assignment[head][0][0] != head_nozzle_assignment[head][-1][0]:
+                    nozzle_change_counter += 1
+            else:
+                if head_nozzle_assignment[head][cycle][0] != head_nozzle_assignment[head][cycle + 1][0]:
+                    nozzle_change_counter += 1
+            head_cycle_counter += head_nozzle_assignment[head][cycle][1]
+        cycle_counter = max(cycle_counter, head_cycle_counter)
+
+    # === 元件拾取次数预估 ===
+    cp_info = []
+    for idx, points in enumerate(assignment_points):
+        if points == 0:
+            continue
+        reminder_points = points % component_feeders[idx]
+        for _ in range(component_feeders[idx]):
+            cp_info.append(
+                [idx, points // component_feeders[idx] + (1 if reminder_points > 0 else 0), component_nozzle[idx]])
+            reminder_points -= 1
+
+    cp_info.sort(key=lambda x: -x[1])
+
+    nozzle_level, nozzle_counter = defaultdict(int), defaultdict(int)
+    level_points = defaultdict(int)
+    for info in cp_info:
+        nozzle = info[2]
+        if nozzle_counter[nozzle] and nozzle_counter[nozzle] % nozzle_heads[nozzle] == 0:
+            nozzle_level[nozzle] += 1
+        level = nozzle_level[nozzle]
+        level_points[level] = max(level_points[level], info[1])
+        nozzle_counter[nozzle] += 1
+
+    pickup_counter = sum(points for points in level_points.values())
+    placement_counter = sum(assignment_points)
+
+    pickup_movement = 0
+    for points in assignment_points:
+        if points:
+            pickup_movement += 1
+    # 返回加权预估时间
+    return t_cycle * cycle_counter + t_nozzle_change * nozzle_change_counter + t_pick * pickup_counter + \
+        t_place * placement_counter + 0.1 * pickup_movement
 
 
-def assemblyline_optimizer_heuristic(pcb_data, component_data):
+def assemblyline_optimizer_heuristic(pcb_data, component_data, machine_number):
     # the number of placement points, the number of available feeders, and nozzle type of component respectively
     component_number = len(component_data)
 
@@ -70,13 +228,15 @@ def assemblyline_optimizer_heuristic(pcb_data, component_data):
         nozzle_points[nozzle] += 1
 
     # first step: generate the initial solution with equalized workload
-    assignment_result = [[0 for _ in range(len(component_points))] for _ in range(max_machine_index)]
-    assignment_points = [0 for _ in range(max_machine_index)]
+    assignment_result = [[0 for _ in range(len(component_points))] for _ in range(machine_number)]
+    assignment_points = [0 for _ in range(machine_number)]
+    average_points = len(pcb_data) // machine_number
 
     weighted_points = list(
         map(lambda x: x[1] + 1e-5 * nozzle_points[component_nozzle[x[0]]], enumerate(component_points)))
 
-    for part_index in np.argsort(weighted_points):
+    # for part_index in np.argsort(weighted_points)[::-1]:
+    for part_index in np.argsort(weighted_points)[::-1]:
         if (total_points := component_points[part_index]) == 0:        # total placements for each component type
             continue
         machine_set = []
@@ -85,7 +245,10 @@ def assemblyline_optimizer_heuristic(pcb_data, component_data):
         for machine_index in np.argsort(assignment_points):
             if len(machine_set) >= component_points[part_index] or len(machine_set) >= component_feeders[part_index]:
                 break
+
             machine_set.append(machine_index)
+            if weighted_points[part_index] + assignment_points[machine_index] < average_points:
+                break
 
         # Allocation of mounting points to available machines according to the principle of equality
         while total_points:
@@ -137,7 +300,141 @@ def assemblyline_optimizer_heuristic(pcb_data, component_data):
                         assignment_result[machine_index][part_index] += 1
                         total_points -= 1
 
+    prev_max_assembly_time, prev_assignment_result = None, None
+    while True:
         # second step: estimate the assembly time for each machine
+        arranged_feeders = defaultdict(list)
+        for machine_index in range(machine_number):
+            arranged_feeders[machine_index] = [0 for _ in range(len(component_data))]
+
+        for part_index in range(len(component_data)):
+            feeder_limit = component_feeders[part_index]    # 总体可用数
+            for machine_index in range(machine_number):
+                if assignment_result[machine_index][part_index] == 0:
+                    continue
+                feeder_limit -= 1
+                # 已分配元件的机器至少安装1把供料器
+                arranged_feeders[machine_index][part_index] = 1
+                assert feeder_limit >= 0
+
+        for part_index in range(len(component_data)):
+            total_feeder_limit = component_feeders[part_index] - sum(
+                [arranged_feeders[machine_index][part_index] for machine_index in range(machine_number)])
+            while total_feeder_limit > 0:
+                max_ratio, max_ratio_machine = None, -1
+                for machine_index in range(machine_number):
+                    if assignment_result[machine_index][part_index] == 0:
+                        continue
+                    ratio = assignment_result[machine_index][part_index] / arranged_feeders[machine_index][part_index]
+                    if max_ratio is None or ratio > max_ratio:
+                        max_ratio, max_ratio_machine = ratio, machine_index
+                assert max_ratio_machine is not None
+                arranged_feeders[max_ratio_machine][part_index] += 1
+                total_feeder_limit -= 1
+
+        assembly_time, chip_per_hour = [], []
+        for machine_index in range(machine_number):
+            assembly_time.append(
+                assembly_time_estimator(assignment_result[machine_index], arranged_feeders[machine_index],
+                                        component_nozzle))
+            chip_per_hour.append(sum(assignment_result[machine_index]) / (assembly_time[-1] + 1e-10))
+
+        max_assembly_time = max(assembly_time)
+        if prev_max_assembly_time and (prev_max_assembly_time < max_assembly_time or abs(
+                max_assembly_time - prev_max_assembly_time) < 1e-10):
+            if prev_max_assembly_time < max_assembly_time:
+                assignment_result = copy.deepcopy(prev_assignment_result)
+            break
+        else:
+            prev_max_assembly_time = max_assembly_time
+            prev_assignment_result = copy.deepcopy(assignment_result)
+
         # third step: adjust the assignment results to reduce maximal assembly time among all machines
+        # ideal averagely assigned points
+        total_points = len(pcb_data)
+        average_assign_points = [round(total_points * chip_per_hour[mi] / sum(chip_per_hour)) for mi in
+                                 range(machine_number)]
+
+        machine_index = 0
+        while total_points != sum(average_assign_points):
+            if total_points > sum(average_assign_points):
+                average_assign_points[machine_index] += 1
+            else:
+                average_assign_points[machine_index] -= 1
+            machine_index += 1
+            if machine_index >= machine_number:
+                machine_index = 0
+
+        # the placement points that need to be re-allocated
+        machine_reallocate_points = [sum(assignment_result[mi]) - average_assign_points[mi] for mi in
+                                     range(machine_number)]
+
+        # workload balance
+        # 1. balance the number of placements of the same type between different machines.
+        for demand_mi in range(machine_number):
+            if machine_reallocate_points[demand_mi] >= 0:
+                continue
+
+            supply_machine_list = [mi for mi in range(machine_number) if machine_reallocate_points[mi] > 0]
+            supply_machine_list.sort(key=lambda mi: -machine_reallocate_points[mi])
+
+            for supply_mi in supply_machine_list:
+                for part_index in range(len(component_data)):
+                    if assignment_result[supply_mi][part_index] <= 0:
+                        continue
+
+                    reallocate_points = min(assignment_result[supply_mi][part_index],
+                                            -machine_reallocate_points[demand_mi])
+
+                    # upper available feeder restrictions
+                    tmp_reallocate_result = [assignment_result[mi][part_index] for mi in range(machine_number)]
+                    tmp_reallocate_result[supply_mi] -= reallocate_points
+                    tmp_reallocate_result[demand_mi] += reallocate_points
+
+                    if sum(1 for pt in tmp_reallocate_result if pt > 0) > component_feeders[part_index]:
+                        continue
+
+                    assignment_result[supply_mi][part_index] -= reallocate_points
+                    machine_reallocate_points[supply_mi] -= reallocate_points
+                    assignment_result[demand_mi][part_index] += reallocate_points
+                    machine_reallocate_points[demand_mi] += reallocate_points
+
+                    if machine_reallocate_points[demand_mi] <= 0:
+                        break
+
+        # 2. balance the number of placements of the different type between different machines.
+        cp_info = []
+        for part_index in range(len(component_data)):
+            for machine_index in range(machine_number):
+                if assignment_result[machine_index][part_index] == 0:
+                    continue
+                cp_info.append([machine_index, part_index, assignment_result[machine_index][part_index]])
+
+        for machine_index in range(machine_number):
+            if machine_reallocate_points[machine_index] >= 0:
+                continue
+
+            filter_cp_info = [info for info in cp_info if
+                              info[0] != machine_index and machine_reallocate_points[info[0]] > 0]
+            while True:
+                if len(filter_cp_info) == 0 or machine_reallocate_points[machine_index] >= 0:
+                    break
+                # todo: 对同时拾取数的影响
+                filter_cp_info.sort(key=lambda x: x[2] + machine_reallocate_points[machine_index])
+
+                info = filter_cp_info[0]
+                filter_cp_info.remove(info)
+                if abs(machine_reallocate_points[machine_index]) + abs(machine_reallocate_points[info[0]]) < abs(
+                        machine_reallocate_points[machine_index] + info[2]) + abs(
+                        machine_reallocate_points[info[0]] - info[2]):
+                    continue
+
+                cp_info.remove(info)
+
+                assignment_result[info[0]][info[1]] = 0
+                assignment_result[machine_index][info[1]] += info[2]
+
+                machine_reallocate_points[info[0]] -= info[2]
+                machine_reallocate_points[machine_index] += info[2]
 
     return assignment_result

optimizer_hyperheuristic.py

@@ -0,0 +1,143 @@
+import pickle
+
+import numpy as np
+
+from base_optimizer.optimizer_interface import *
+from generator import *
+
+os.environ['KMP_DUPLICATE_LIB_OK'] = 'True'
+
+
+class Net(torch.nn.Module):
+    def __init__(self, input_size, output_size):
+        super(Net, self).__init__()
+        self.fc1 = torch.nn.Linear(input_size, 1024)
+        self.relu = torch.nn.ReLU()  # 激活函数
+        self.fc2 = torch.nn.Linear(1024, output_size)
+
+    def forward(self, x):
+        x = self.fc1(x)
+        x = self.relu(x)
+        x = self.fc2(x)
+        return x
+
+
+def selective_initialization(component_points, population_size, machine_number):
+    # assignment_result = [[0 for _ in range(len(component_points))] for _ in range(machine_number)]
+    assignment_result = []
+
+    return assignment_result
+
+
+def optimizer_hyperheuristc(pcb_data, component_data, machine_number):
+
+    # genetic-based hyper-heuristic
+    crossover_rate, mutation_rate = 0.8, 0.1
+    population_size, n_generations = 200, 500
+
+    # todo: how to generate initial population (random?)
+    # assignment_result = selective_initialization(component_points, population_size, machine_number)
+    assignment_result = []
+    return assignment_result
+
+
+if __name__ == '__main__':
+    warnings.simplefilter(action='ignore', category=FutureWarning)
+
+    train_file, test_file = 'train_data.txt', 'test_data.txt'
+    num_epochs = 30000
+
+    data_mgr = DataMgr()
+    device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
+
+    # batch_size = 40000
+    # for _ in range(batch_size):
+    #     pcb_data, component_data = data_mgr.generator()     # random generate a PCB data
+    #     # data_mgr.remover()                                # 移除最近一次保存数据
+    #     # data_mgr.saver('data/' + train_file, pcb_data)     # 保存新数据
+    #
+    #     info = base_optimizer(1, pcb_data, component_data,
+    #                           feeder_data=pd.DataFrame(columns=['slot', 'part', 'arg']), method='feeder_scan',
+    #                           hinter=True)
+    #
+    #     data_mgr.recorder('opt/' + train_file, info, pcb_data, component_data)
+
+    train, save = True, True
+    learning_rate = 0.0005
+
+    net = Net(input_size=data_mgr.get_feature(), output_size=1).to(device)
+    if train:
+        data = data_mgr.loader('opt/' + train_file)
+        x_fit, y_fit = np.array(data[2:]).T, np.array([data[1]]).T
+        lr = LinearRegression()
+        lr.fit(x_fit, y_fit)
+
+        x_train, y_train = np.array(data[0]), lr.predict(x_fit)     # np.array(data[1])
+        x_train = torch.from_numpy(x_train.reshape((-1, np.shape(x_train)[1]))).float().to(device)
+        y_train = torch.from_numpy(y_train.reshape((-1, 1))).float().to(device)
+
+        optimizer = torch.optim.Adam(net.parameters(), lr=learning_rate)
+        # scheduler = torch.optim.lr_scheduler.StepLR(optimizer, step_size=1000, gamma=0.1)
+
+        loss_func = torch.nn.MSELoss()
+
+        for epoch in range(num_epochs):
+            pred = net(x_train)
+            loss = loss_func(pred, y_train)
+            optimizer.zero_grad()
+            loss.backward()
+            optimizer.step()
+            # scheduler.step()
+            if epoch % 200 == 0:
+                print('Epoch:  ', epoch, ', Loss: ', loss.item())
+                if loss.item() < 1e-4:
+                    break
+
+        net_predict = net(x_train).view(-1)
+        # pred_time, real_time = net_predict.cpu().detach().numpy(), y_train.view(-1).cpu().detach().numpy()
+        pred_time, real_time = net_predict.cpu().detach().numpy(), np.array(data[1])
+
+        pred_error = np.array([])
+        for t1, t2 in np.nditer([pred_time, real_time]):
+            pred_error = np.append(pred_error, abs(t1 - t2) / (t2 + 1e-10) * 100)
+
+        print('--------------------------------------')
+        print(f'average prediction error for train data : {np.average(pred_error): .2f}% ')
+        print(f'maximum prediction error for train data : {np.max(pred_error): .2f}% ')
+
+        mse = np.linalg.norm((net_predict - y_train.view(-1)).cpu().detach().numpy())
+        print(f'mean square error for training data result : {mse: 2f} ')
+        if save:
+            torch.save(net.state_dict(), 'model_state.pth')
+            with open('lr_model.pkl', 'wb') as f:
+                pickle.dump(lr, f)
+            joblib.dump(lr, "lr_model.m")
+            # torch.save(optimizer.state_dict(), 'optimizer_state.pth')
+    else:
+        with open('lr_model.pkl', 'rb') as f:
+            lr = pickle.load(f)
+        net.load_state_dict(torch.load('model_state.pth'))
+        # optimizer = torch.optim.Adam(net.parameters(), lr=learning_rate)
+        # optimizer.load_state_dict(torch.load('optimizer_state.pth'))
+
+    data = data_mgr.loader('opt/' + test_file)
+    x_test, y_test = np.array(data[0]), lr.predict(np.array(data[2:]).T)
+    x_test, y_test = torch.from_numpy(x_test.reshape((-1, np.shape(x_test)[1]))).float().to(device), \
+                     torch.from_numpy(y_test.reshape((-1, 1))).float().to(device)
+
+    net.eval()
+    with torch.no_grad():
+        net_predict = net(x_test).view(-1)
+        pred_time, real_time = net_predict.cpu().detach().numpy(), y_test.view(-1).cpu().detach().numpy()
+
+        pred_error = np.array([])
+        for t1, t2 in np.nditer([pred_time, real_time]):
+            pred_error = np.append(pred_error, abs(t1 - t2) / (t2 + 1e-10) * 100)
+
+        print('--------------------------------------')
+        print(f'average prediction error for test data : {np.average(pred_error): .2f}% ')
+        print(f'maximum prediction error for test data : {np.max(pred_error): .2f}% ')
+
+        mse = np.linalg.norm(pred_time - real_time)
+        print(f'mean square error for test data result : {mse: 2f} ')
+

optimizer_reconfiguration.py

@@ -1,12 +1,8 @@
-import copy
-import math
-import random
-import numpy as np
-
 from base_optimizer.optimizer_common import *
 
 
-def objective_value_calculate(component_assignment, component_nozzle, task_block_weight):
+# 生产过程中不允许吸嘴更换/点的拾取贴装仅与供料器槽位/模组相关
+def objective_value_calculate(component_assignment, component_nozzle, task_block_weight, machine_number):
     machine_assembly_time = []
     for machine_index in range(max_machine_index):
         task_block_number, total_point_number = 0, sum(component_assignment[machine_index])
@@ -39,10 +35,10 @@ def objective_value_calculate(component_assignment, component_nozzle, task_block
     return max(machine_assembly_time)
 
 
-def random_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight):
+def random_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight, machine_number):
     component_points_cpy = copy.deepcopy(component_points)
     component_number = len(component_points_cpy)
-    assignment_result = [[0 for _ in range(component_number)] for _ in range(max_machine_index)]
+    assignment_result = [[0 for _ in range(component_number)] for _ in range(machine_number)]
 
     # == the set of feasible component type for each nozzle type
     nozzle_part_list = defaultdict(list)
@@ -52,7 +48,7 @@ def random_component_assignment(component_points, component_nozzle, component_fe
     selected_part = []
     for part_list in nozzle_part_list.values():
         part = random.sample(part_list, 1)[0]
-        machine_index = random.randint(0, max_machine_index - 1)
+        machine_index = random.randint(0, machine_number - 1)
 
         assignment_result[machine_index][part] += 1
         component_points_cpy[part] -= 1
@@ -63,24 +59,24 @@ def random_component_assignment(component_points, component_nozzle, component_fe
         if part in selected_part:
             continue
 
-        assignment_result[random.randint(0, max_machine_index - 1)][part] += 1
+        assignment_result[random.randint(0, machine_number - 1)][part] += 1
         component_points_cpy[part] -= 1
 
-    machine_assign = list(range(max_machine_index))
+    machine_assign = list(range(machine_number))
     random.shuffle(machine_assign)
     finished_assign_counter = 0
     while finished_assign_counter < component_number:
-        # todo: feeder limit restriction
         for machine_index in machine_assign:
             part = random.randint(0, component_number - 1)
             feeder_counter = 0
-            for idx in range(max_machine_index):
+            for idx in range(machine_number):
                 if assignment_result[idx][part] > 0 or idx == machine_index:
                     feeder_counter += 1
 
             if component_points_cpy[part] == 0 or feeder_counter > component_feeders[part]:
                 continue
 
+            # feeder limit restriction
             points = random.randint(1, component_points_cpy[part])
             assignment_result[machine_index][part] += points
             component_points_cpy[part] -= points
@@ -96,15 +92,16 @@ def greedy_component_assignment(component_points, component_nozzle, component_fe
     pass    # 不清楚原文想说什么
 
 
-def local_search_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight):
+def local_search_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight,
+                                      machine_number):
     # maximum number of iterations : 5000
     # maximum number of unsuccessful iterations: 50
     component_number = len(component_points)
     iteration_counter, unsuccessful_iteration_counter = 5000, 50
     optimal_val, optimal_assignment = random_component_assignment(component_points, component_nozzle, component_feeders,
-                                                                  task_block_weight)
+                                                                  task_block_weight, machine_number)
     for _ in range(iteration_counter):
-        machine_index = random.randint(0, max_machine_index - 1)
+        machine_index = random.randint(0, machine_number - 1)
         if sum(optimal_assignment[machine_index]) == 0:
             continue
 
@@ -120,9 +117,9 @@ def local_search_component_assignment(component_points, component_nozzle, compon
         swap_machine_index = None
         while cyclic_counter <= 2 * machine_index:
             cyclic_counter += 1
-            swap_machine_index = random.randint(0, max_machine_index - 1)
+            swap_machine_index = random.randint(0, machine_number - 1)
             feeder_available = 0
-            for machine in range(max_machine_index):
+            for machine in range(machine_number):
                 if optimal_assignment[machine][component_index] or machine == swap_machine_index:
                     feeder_available += 1
 
@@ -143,18 +140,18 @@ def local_search_component_assignment(component_points, component_nozzle, compon
     return optimal_val, optimal_assignment
 
 
-def reconfig_crossover_operation(component_points, component_feeders, parent1, parent2):
+def reconfig_crossover_operation(component_points, component_feeders, parent1, parent2, machine_number):
     offspring1, offspring2 = copy.deepcopy(parent1), copy.deepcopy(parent2)
     component_number = len(component_points)
 
     # === crossover ===
     mask_bit = []
-    for _ in range(max_machine_index):
+    for _ in range(machine_number):
         mask_bit.append(random.randint(0, 1))
-    if sum(mask_bit) == 0 or sum(mask_bit) == max_machine_index:
+    if sum(mask_bit) == 0 or sum(mask_bit) == machine_number:
         return offspring1, offspring2
 
-    for machine_index in range(max_machine_index):
+    for machine_index in range(machine_number):
         if mask_bit:
             offspring1[machine_index] = copy.deepcopy(parent1[machine_index])
             offspring2[machine_index] = copy.deepcopy(parent2[machine_index])
@@ -166,20 +163,20 @@ def reconfig_crossover_operation(component_points, component_feeders, parent1, p
     # equally to reach the correct number
     for component_index in range(component_number):
         for offspring in [offspring1, offspring2]:
-            additional_points = sum([offspring[mt][component_index] for mt in range(max_machine_index)]) - \
+            additional_points = sum([offspring[mt][component_index] for mt in range(machine_number)]) - \
                                 component_points[component_index]
             if additional_points > 0:
                 # if a component type has more placements, decrease the assigned values on every head equally keeping
                 # the proportion of the number of placement among the heads
                 points_list = []
-                for machine_index in range(max_machine_index):
+                for machine_index in range(machine_number):
                     points = math.floor(
                         additional_points * offspring[machine_index][component_index] / component_points[component_index])
                     points_list.append(points)
                     offspring[machine_index][component_index] -= points
                 additional_points -= sum(points_list)
 
-                for machine_index in range(max_machine_index):
+                for machine_index in range(machine_number):
                     if additional_points == 0:
                         break
                     if offspring[machine_index][component_index] == 0:
@@ -189,7 +186,7 @@ def reconfig_crossover_operation(component_points, component_feeders, parent1, p
             elif additional_points < 0:
                 # otherwise, increase the assigned nonzero values equally
                 machine_set = []
-                for machine_index in range(max_machine_index):
+                for machine_index in range(machine_number):
                     if offspring[machine_index][component_index] == 0:
                         continue
                     machine_set.append(machine_index)
@@ -205,29 +202,23 @@ def reconfig_crossover_operation(component_points, component_feeders, parent1, p
                     offspring[machine_index][component_index] += 1
                     additional_points -= 1
 
+    # === 结果校验 ===
+    for offspring in [offspring1, offspring2]:
         for part in range(component_number):
-        pt = 0
-        for mt in range(max_machine_index):
-            pt+= offspring1[mt][part]
-        if pt!=component_points[part]:
-            print('')
-    for part in range(component_number):
-        pt = 0
-        for mt in range(max_machine_index):
-            pt+= offspring2[mt][part]
-        if pt!=component_points[part]:
-            print('')
+            pt = sum(offspring[mt][part] for mt in range(machine_number))
+            assert pt == component_points[part]
+
     return offspring1, offspring2
 
 
-def reconfig_mutation_operation(component_feeders, parent):
+def reconfig_mutation_operation(component_feeders, parent, machine_number):
     offspring = copy.deepcopy(parent)
 
     swap_direction = random.randint(0, 1)
     if swap_direction:
-        swap_machine1, swap_machine2 = random.sample(list(range(max_machine_index)), 2)
+        swap_machine1, swap_machine2 = random.sample(list(range(machine_number)), 2)
     else:
-        swap_machine2, swap_machine1 = random.sample(list(range(max_machine_index)), 2)
+        swap_machine2, swap_machine1 = random.sample(list(range(machine_number)), 2)
 
     component_list = []
     for component_index, points in enumerate(offspring[swap_machine1]):
@@ -248,7 +239,7 @@ def reconfig_mutation_operation(component_feeders, parent):
     return offspring
 
 
-def evolutionary_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight):
+def evolutionary_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight, machine_number):
     # population size: 10
     # probability of the mutation: 0.1
     # probability of the crossover: 0.8
@@ -260,7 +251,8 @@ def evolutionary_component_assignment(component_points, component_nozzle, compon
     population = []
     for _ in range(population_size):
         population.append(
-            random_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight)[1])
+            random_component_assignment(component_points, component_nozzle, component_feeders, task_block_weight,
+                                        machine_number)[1])
 
     with tqdm(total=generation_number) as pbar:
         pbar.set_description('evolutionary algorithm process for PCB assembly line balance')
@@ -270,7 +262,7 @@ def evolutionary_component_assignment(component_points, component_nozzle, compon
             # calculate fitness value
             pop_val = []
             for individual in population:
-                pop_val.append(objective_value_calculate(individual, component_nozzle, task_block_weight))
+                pop_val.append(objective_value_calculate(individual, component_nozzle, task_block_weight, machine_number))
 
             select_index = get_top_k_value(pop_val, population_size - len(new_population), reverse=False)
             population = [population[idx] for idx in select_index]
@@ -278,7 +270,7 @@ def evolutionary_component_assignment(component_points, component_nozzle, compon
 
             population += new_population
             for individual in new_population:
-                pop_val.append(objective_value_calculate(individual, component_nozzle, task_block_weight))
+                pop_val.append(objective_value_calculate(individual, component_nozzle, task_block_weight, machine_number))
 
             # min-max convert
             max_val = max(pop_val)
@@ -297,13 +289,14 @@ def evolutionary_component_assignment(component_points, component_nozzle, compon
                             break
 
                     offspring1, offspring2 = reconfig_crossover_operation(component_points, component_feeders,
-                                                                          population[index1], population[index2])
+                                                                          population[index1], population[index2],
+                                                                          machine_number)
 
                     if np.random.random() < mutation_rate:
-                        offspring1 = reconfig_mutation_operation(component_feeders, offspring1)
+                        offspring1 = reconfig_mutation_operation(component_feeders, offspring1, machine_number)
 
                     if np.random.random() < mutation_rate:
-                        offspring2 = reconfig_mutation_operation(component_feeders, offspring2)
+                        offspring2 = reconfig_mutation_operation(component_feeders, offspring2, machine_number)
 
                     new_population.append(offspring1)
                     new_population.append(offspring2)
@@ -313,7 +306,7 @@ def evolutionary_component_assignment(component_points, component_nozzle, compon
     return min(pop_val), population[np.argmin(pop_val)]
 
 
-def reconfiguration_optimizer(pcb_data, component_data):
+def reconfiguration_optimizer(pcb_data, component_data, machine_number):
     # === data preparation ===
     component_number = len(component_data)
 
@@ -335,20 +328,28 @@ def reconfiguration_optimizer(pcb_data, component_data):
     # === assignment of heads to modules is omitted ===
     optimal_assignment, optimal_val = [], None
 
-    task_block_weight = 5   # element from list [0, 1, 2, 5, 10]
+    task_block_weight = 5   # element from list [0, 1, 2, 5, 10] task_block ~= cycle
     # === assignment of components to heads
-    for i in range(4):
+    for i in range(5):
         if i == 0:
+            # random
             val, assignment = random_component_assignment(component_points, component_nozzle, component_feeders,
-                                                          task_block_weight)
+                                                          task_block_weight, machine_number)
         elif i == 1:
+            # brute force
+            # which is proved to be useless, since it only ran in reasonable time for the smaller test instances
             continue
         elif i == 2:
-            val, assignment = local_search_component_assignment(component_points, component_nozzle,
-                                                                component_feeders, task_block_weight)
+            # local search
+            val, assignment = local_search_component_assignment(component_points, component_nozzle, component_feeders,
+                                                                task_block_weight, machine_number)
+        elif i == 3:
+            # evolutionary
+            val, assignment = evolutionary_component_assignment(component_points, component_nozzle, component_feeders,
+                                                                task_block_weight, machine_number)
         else:
-            val, assignment = evolutionary_component_assignment(component_points, component_nozzle,
-                                                                component_feeders, task_block_weight)
+            # greedy: unclear description
+            continue
 
         if optimal_val is None or val < optimal_val:
             optimal_val, optimal_assignment = val, assignment.copy()