OpenMP 旅行商问题,静态调度

▶ 《并行程序设计导论》第六章中讨论了旅行商,分别使用了 MPI,Pthreads,OpenMP 来进行实现,这里是 OpenMP 的代码,分为静态调度(每个线程分分配等量的搜索人物)和动态调度(每个线程分配不等量的任务,每当有线程完成自己的任务后,向其他线程请求新的子任务)

● 静态调度代码

 #include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <omp.h> #define DEBUG
#define City_count(tour) (tour->count)
#define Tour_cost(tour) (tour->cost)
#define Last_city(tour) (tour->cityList[(tour->count)-1]) // tour 中最后一个城市的编号
#define Tour_city(tour,i) (tour->cityList[(i)]) // tour 中第 i 个城市的编号
#define Cost(city1, city2) (digraph[city1 * n + city2]) // tour 中两个城市之间的代价
#define Queue_elt(queue,i) (queue->list[(queue->head + (i)) % queue->spaceAlloc]) typedef int city_t; // 城市编号
typedef int cost_t; // 代价,可以为浮点数 typedef struct
{
int count; // 已经走过的城市数
city_t* cityList; // 已经走过的城市列表
cost_t cost; // 当前代价
} tour_struct; typedef struct
{
int spaceAlloc; // 栈占据的空间大小(最大大小)
int stackSize; // 栈的实际大小
tour_struct** list;
} stack_struct; typedef struct
{
int spaceAlloc; // 队列占据的空间大小(最大大小)
int head; // 队头
int tail; // 队尾
bool full; // 是否队满
tour_struct** list;
} queue_struct; const int INFINITY = ; // 最大代价,用于初始化 tour
const int NO_CITY = -; // 标记 tour->cityList 中尚未使用的部分
const int MAX_STRING = ; // 输出字符串最大长度
int n; // 城市总数
int nThread; // 线程数
city_t home_town = ; // 起点和终点的编号
cost_t *digraph; // 读取的邻接表数据
tour_struct *best_tour; // 保存最优轨迹
queue_struct *queue; // 搜索队列
int queueSize; // 搜索队的大小
int initialTourCount; // 首次分配给各线程的搜索队列数量
tour_struct* allocTour(stack_struct* avail);// 必需的声明,不然该函数与栈操作函数混在一起 // 基本函数
void usage(char* prog_name)// 提示信息
{
fprintf(stderr, "\n<Usage> %s <nThread> <digraph file>\n", prog_name);
exit();
} void readDigraph(FILE* digraph_file)// 从文件读取邻接表
{
int i, j;
fscanf_s(digraph_file, "%d", &n);// 第一行数字时城市数量
if (n <= )
{
fprintf(stderr, "\n<readDigraph> Number of city = %d\n", n);
exit();
}
digraph = (cost_t*)malloc(n * n * sizeof(cost_t));
for (i = ; i < n; i++)
{
for (j = ; j < n; j++)
{
fscanf_s(digraph_file, "%d", &digraph[i * n + j]);
if (i == j && digraph[i * n + j] != )
{
fprintf(stderr, "\n<readDigraph> Diagonal element [%d, %d] = %d\n", i, j, digraph[i * n + j]);
exit();
}
else if (i != j && digraph[i * n + j] <= )
{
fprintf(stderr, "\n<readDigraph> Off-diagonal element [%d, %d] = %d\n", i, j, digraph[i * n + j]);
exit();
}
}
}
#ifdef DEBUG
printf("\n\t<readDigraph> Read map\n");
#endif
} void printDigraph(void)// 打印邻接表
{
int i, j;
printf("\n\t<printDigraph> Order = %d\nMatrix = \n", n);
for (i = ; i < n; i++)
{
for (j = ; j < n; j++)
printf("%2d ", digraph[i * n + j]);
printf("\n");
}
printf("\n");
} void initeTour(tour_struct* tour, cost_t cost)// 初始化搜索轨迹,要求给出代价
{
int i;
tour->cityList[] = ;
for (i = ; i <= n; tour->cityList[i++] = NO_CITY);
tour->cost = cost;
tour->count = ;
} void copyTour(tour_struct* tour1, tour_struct* tour2)// 搜索轨迹 tour1 拷贝给 tour2
{
memcpy(tour2->cityList, tour1->cityList, (n + ) * sizeof(city_t));
tour2->count = tour1->count;
tour2->cost = tour1->cost;
} void printTour(int my_rank, tour_struct* tour, char* title)// 输出搜索轨迹,注意有两种输出
{
int i;
char string[MAX_STRING];
if (my_rank >= )// 各线程汇报
sprintf_s(string, "\n\t<printTour> Rank%2d, %s %p: ", my_rank, title, tour);
else // 输出当前最佳回路
sprintf_s(string, "\n\t<printTour> %s: ", title);
for (i = ; i < City_count(tour); i++)
sprintf_s(string + strlen(string), , "%d ", Tour_city(tour, i));
printf("%s\n", string);
} // 栈相关
stack_struct* initeStack(void)// 初始化栈
{
int i;
stack_struct* stack = (stack_struct*)malloc(sizeof(stack_struct));
stack->list = (tour_struct**)malloc(n * n * sizeof(tour_struct*));
for (i = ; i < n * n; stack->list[i++] = NULL);
stack->stackSize = ;
stack->spaceAlloc = n * n;
return stack;
} void pushStack(stack_struct* stack, tour_struct* tour)// 非拷贝压栈
{
if (stack->stackSize == stack->spaceAlloc)
{
fprintf(stderr, "\n<pushStack> Overflow\n");
free(tour->cityList);
free(tour);
}
else
{
# ifdef DEBUG
printf("\n\t<pushStack> StackSize = %d, tour at %p, tour->cityList at %p\n", stack->stackSize, tour, tour->cityList);
printTour(-, tour, "pushStack");
printf("\n");
# endif
stack->list[stack->stackSize] = tour;
(stack->stackSize)++;
}
} void pushCopyStack(stack_struct* stack, tour_struct* tour, stack_struct* avail)// 拷贝压栈
{
tour_struct* tmp;
if (stack->stackSize == stack->spaceAlloc)
{
fprintf(stderr, "\n<pushCopyStack> Overflow\n");
exit(-);
}
tmp = allocTour(avail);
copyTour(tour, tmp);
stack->list[stack->stackSize] = tmp;
(stack->stackSize)++;
} tour_struct* popStack(stack_struct* stack)// 出栈
{
tour_struct* tmp;
if (stack->stackSize == )
{
fprintf(stderr, "\n<popStack> Empty\n");
exit(-);
}
tmp = stack->list[stack->stackSize - ];
stack->list[stack->stackSize - ] = NULL;
(stack->stackSize)--;
return tmp;
} int isEmptyStack(stack_struct* stack)// 判断是否空栈
{
return stack->stackSize == ;
} void freeStack(stack_struct* stack)// 清空栈
{
int i;
for (i = ; i < stack->stackSize; free(stack->list[i]->cityList), free(stack->list[i]), i++);
free(stack->list);
free(stack);
} void printStack(stack_struct* stack, int my_rank, char title[])// 打印栈
{
char string[MAX_STRING];
int i, j;
printf("\n\t<printStack> Rank%2d, %s\n", my_rank, title);
for (i = ; i < stack->stackSize; i++)
{
sprintf_s(string, , "%d> ", i);
for (j = ; j < stack->list[i]->count; j++)
sprintf_s(string + strlen(string), , "%d ", stack->list[i]->cityList[j]);
printf("%s\n", string);
}
} // 队列相关
queue_struct* initeQueue(int size)// 初始化队列
{
queue_struct* new_queue = (queue_struct*)malloc(sizeof(queue_struct));
new_queue->list = (tour_struct**)malloc(size * sizeof(tour_struct*));
new_queue->spaceAlloc = size;
new_queue->head = new_queue->tail = new_queue->full = false;
return new_queue;
} int isEmptyQueue(queue_struct* queue)// 判断是否队空
{
return !queue->full && queue->head == queue->tail;// 空队要求 queue->full == false 且头尾指针相等
} tour_struct* deQueue(queue_struct* queue)// 出队
{
tour_struct* tmp;
if (isEmptyQueue(queue))
{
fprintf(stderr, "\n<deQueue> Empty queue\n");
exit(-);
}
tmp = queue->list[queue->head];
queue->head = (queue->head + ) % queue->spaceAlloc;
return tmp;
} void enQueue(queue_struct* queue, tour_struct* tour)// 入队
{
tour_struct* tmp;
if (queue->full == true)
{
fprintf(stderr, "\n<enQueue> Overflow\n");
exit(-);
}
tmp = allocTour(NULL);
copyTour(tour, tmp);
queue->list[queue->tail] = tmp;
queue->tail = (queue->tail + ) % queue->spaceAlloc;
if (queue->tail == queue->head)
queue->full = true;
} void freeQueue(queue_struct* queue)// 清空队列
{
free(queue->list);
free(queue);
} void printQueue(queue_struct* queue, int my_rank, char title[])// 打印队列
{
char string[MAX_STRING];
int i, j;
printf("\n\t<printQueue> Rank%2d > %s\n", my_rank, title);
for (i = queue->head; i != queue->tail; i = (i + ) % queue->spaceAlloc)
{
sprintf_s(string, "%d> %p = ", i, queue->list[i]);
for (j = ; j < queue->list[i]->count; j++)
sprintf_s(string + strlen(string), , "%d ", queue->list[i]->cityList[j]);
printf("%s\n", string);
}
} // Tour 相关
tour_struct* allocTour(stack_struct* avail)// 生成一个搜索轨迹
{
tour_struct* tmp;
if (avail == NULL || isEmptyStack(avail))
{
tmp = (tour_struct*)malloc(sizeof(tour_struct));
tmp->cityList = (city_t*)malloc((n + ) * sizeof(city_t));
return tmp;
}
return popStack(avail);
} int visited(tour_struct* tour, city_t city)// 判断一个轨迹中是否经过了城市 city
{
for (int i = ; i < City_count(tour); i++)
{
if (Tour_city(tour, i) == city)
return true;
}
return false;
} void addCity(tour_struct* tour, city_t new_city)// 将一个城市添加到 tour
{
city_t old_last_city = Last_city(tour); // 从原来的最后一个城市接出一条路径
tour->cityList[tour->count] = new_city;
(tour->count)++;
tour->cost += Cost(old_last_city, new_city);
} void removeLastCity(tour_struct* tour)// 去掉 tour 中最后一个城市
{
city_t old_last_city = Last_city(tour), new_last_city;
tour->cityList[tour->count - ] = NO_CITY;
(tour->count)--;
new_last_city = Last_city(tour);
tour->cost -= Cost(new_last_city, old_last_city);
} int isBestTour(tour_struct* tour)// 判断当前轨迹是否最佳
{
cost_t cost_so_far = Tour_cost(tour);
city_t last_city = Last_city(tour);
return cost_so_far + Cost(last_city, home_town) < Tour_cost(best_tour);// 当前轨迹的代价加上最后一个城市回家的代价是否更优
} int isFeasible(tour_struct* tour, city_t city)// 判断 tour 的最后一个城市是否可以前往城市 city
{
return !visited(tour, city) && Tour_cost(tour) + Cost(Last_city(tour), city) < Tour_cost(best_tour);// 要求没去过该城市且代价不超出当前最优回路代价
} void updateBestTour(tour_struct* tour)
{
if (isBestTour(tour))// 再次检查是否最优,防止两次检查之间其他线程更新了最优轨迹
{
copyTour(tour, best_tour);
addCity(best_tour, home_town);
}
} void freeTour(tour_struct* tour, stack_struct* avail)
{
if (avail == NULL) // 普通释放
{
free(tour->cityList);
free(tour);
}
else // 将结点压回栈中
pushStack(avail, tour);
} void distributeTour(int my_rank, int* my_first_tour_p, int* my_last_tour_p)// 分配初始搜索任务,使用块划分
{
const int quotient = initialTourCount / nThread, remainder = initialTourCount % nThread;// 平均每线程 quotient 个,余数部分由靠前的线程分担
int my_count;
if (my_rank < remainder)// 靠前的线程,分担余数
{
my_count = quotient + ;
*my_first_tour_p = my_rank*my_count;
}
else // 靠后的线程,直接分配,注意偏移
{
my_count = quotient;
*my_first_tour_p = my_rank*my_count + remainder;
}
*my_last_tour_p = *my_first_tour_p + my_count - ;
} // 核心计算函数附属
inline long long factorial(int k)// 计算 k 的阶乘
{
long long tmp = ;
int i;
for (i = ; i <= k; tmp *= i++);
return tmp;
} inline int supQueueSize(void)// 估计搜索需要的队列长度上限,并判断使用的线程数是否太多,没看懂
{
int fact, size;
for (fact = size = n - ; size < nThread; size *= ++fact);
if (size > factorial(n - ))
{
fprintf(stderr, "\n<supQueueSize> Too many thread\n");
size = ;
}
return size;
} void buildInitialQueue(void)// 生成初始搜索队列
{
int currentQueueSize; // 当前队列的长度
city_t nbr;
tour_struct* tour = allocTour(NULL); initeTour(tour, ); // 仅包含起点,代价为 0
queue = initeQueue( * queueSize); // 搜索队 enQueue(queue, tour); // 将起点拷贝入队
freeTour(tour, NULL);
for (currentQueueSize = ; queueSize < nThread;)// 向搜索队中添加元素,直到不小于线程数为止,以便分配各线程工作
{
tour = deQueue(queue); // 从队列中取出一个城市来
currentQueueSize--;
for (nbr = ; nbr < n; nbr++)
{
if (!visited(tour, nbr)) // 取出的城市是新城市
{
addCity(tour, nbr); // 将该城市添加进 tour
enQueue(queue, tour); // 将当前 tour 拷贝加入搜索队列中
currentQueueSize++;
removeLastCity(tour); // 加入搜索队之后立即从 tour 中清楚刚加入的城市,tour 仅用于搜索队操作
}
}
freeTour(tour, NULL); // 本轮入队结束
}
initialTourCount = currentQueueSize;// 记录当前搜索队中 tour 的数量
# ifdef DEBUG
printQueue(queue, , "\n\t<buildInitialQueue> Queue Initialized\n");
# endif
} void treePartition(int my_rank, stack_struct* stack)// 分树
{
int my_first_tour, my_last_tour, i;
# pragma omp single // 估计搜索队列长度上限
queueSize = supQueueSize();
# ifdef DEBUG
printf("\n\t<treePartition> Rank%2d> supQueueSize = %d\n", my_rank, queueSize);
# endif
if (queueSize == ) // 线程数太多
exit(); # pragma omp master // 主线程创建初始搜索队列
buildInitialQueue();
# pragma omp barrier // 从线程要显式的等待主线程 distributeTour(my_rank, &my_first_tour, &my_last_tour);// 分配初始搜索任务
# ifdef DEBUG
printf("\n\t<treePartition> Rank%2d> initialTourCount = %d, first = %d, last = %d\n", my_rank, initialTourCount, my_first_tour, my_last_tour);
# endif
for (i = my_last_tour; i >= my_first_tour; i--)// 从分配到的任务列中中最后一个向前逐个压入栈中
{
# ifdef DEBUG
printTour(my_rank, Queue_elt(queue, i), "pushStack");
# endif
pushStack(stack, Queue_elt(queue, i));
}
# ifdef DEBUG
printStack(stack, my_rank, "Task set up");
# endif } // 核心计算函数
void treeSearch(void)
{
int my_rank = omp_get_thread_num();
city_t nbr;
stack_struct *stack = initeStack(), *avail = initeStack();// 搜索栈和未搜索的栈
tour_struct *curr_tour; for (treePartition(my_rank, stack); !isEmptyStack(stack);)// 初次分树,然后全部搜索完以前不断循环
{
curr_tour = popStack(stack);
# ifdef DEBUG
printTour(my_rank, curr_tour, "popStack");
# endif
if (City_count(curr_tour) == n) // 当前回路已经搜索了 n 个城市
{
if (isBestTour(curr_tour)) // 判断是否是最佳回路
{
# ifdef DEBUG
printTour(my_rank, curr_tour, "Best tour");
# endif
# pragma omp critical
updateBestTour(curr_tour);// 更新当前最佳回路
}
}
else // 还没有搜索 n 个城市
{
for (nbr = n - ; nbr >= ; nbr--)
{
if (isFeasible(curr_tour, nbr))
{
addCity(curr_tour, nbr);
pushCopyStack(stack, curr_tour, avail);
removeLastCity(curr_tour);
}
}
}
freeTour(curr_tour, avail);
}
freeStack(stack); // 释放线程栈资源
freeStack(avail);
# pragma omp barrier // 等待所有线程都完成
# pragma omp master
freeQueue(queue); // 主线程释放队列资源
} int main(int argc, char* argv[])
{
FILE* digraph_file; // 输入文件名
double start, finish; // 计时器
# ifdef DEBUG
argc = ;
argv[] = "";
argv[] = "D:\\Code\\并行程序设计导论 - 代码\\ch6\\TSP\\mat_17e";
# endif if (argc != ) // 检查输入
usage(argv[]);
if ((nThread = strtol(argv[], NULL, )) <= )
{
fprintf(stderr, "\n<main> Error thread number\n");
usage(argv[]);
}
fopen_s(&digraph_file, argv[], "r");
if (digraph_file == NULL)
{
fprintf(stderr, "\n<main> Error opening input file\n");
usage(argv[]);
} readDigraph(digraph_file);// 读取图
fclose(digraph_file);
# ifdef DEBUG
printDigraph();
# endif best_tour = allocTour(NULL);
initeTour(best_tour, INFINITY);
# ifdef DEBUG
printTour(-, best_tour, "Best tour");
printf("\n\t<main> City count = %d\nCost = %d\n\n", City_count(best_tour), Tour_cost(best_tour));
# endif start = omp_get_wtime();// 开始计算
# pragma omp parallel num_threads(nThread) default(none)
treeSearch();
finish = omp_get_wtime(); printTour(-, best_tour, "Best tour");// 汇报结果
printf("\n\t<main> Cost = %d\nElapsed time = %e s\n", best_tour->cost, finish - start); free(best_tour->cityList);// 释放资源
free(best_tour);
free(digraph);
return ;
}

● 输入的邻接表图


● 输出结果,因为是指数级的枚举,没有等待程序运行结束

上一篇:eclipse mars 4.5.1 自定义工具栏


下一篇:Django自定义过滤器中is_safe和need_autoescape两个参数的理解