在 C++ 中通过网格/矩阵找到成本优化路径

Posted

技术标签:

【中文标题】在 C++ 中通过网格/矩阵找到成本优化路径【英文标题】:find a cost optimized path through a grid/matrix in c++ 【发布时间】:2011-09-12 17:17:11 【问题描述】:

我遇到了一个问题,在网上找不到太多帮助。我需要从多个数字向量中找到数字的最小成本组合。所有向量的向量大小相同。 例如,考虑以下内容:

row [0]:  a  b  c  d   
row [1]:  e  f  g  h  
row [2]:  i  j  k  l  

现在我需要通过从每一行(即向量)中取一个元素来找到数字的组合,例如:aei 在此之后,我需要找到其他三个不相交的组合,例如:bfj、cgk、dhl。我根据选择的这四种组合计算成本。目标是找到成本最低的组合。另一种可能的组合可以是:afj、bei、chk、dgl。如果总列数为 d,行数为 k,则可能的总组合为 d^k。这些行存储为向量。我被困在这里,我发现很难为上述过程编写算法。如果有人可以提供帮助,我将不胜感激。 谢谢。

// I am still working on the algorithm. I just have the vectors and the cost function.  

//Cost Function  , it also depends on the path chosen
float cost(int a, int b, PATH to_a)   
float costValue;  
...  
...  
return costValue;  
  

vector< vector < int > > row;  
//populate row  
...   
...
//Suppose  

//    row [0]:  a  b  c  d   
//    row [1]:  e  f  g  h  
//    row [2]:  i  j  k  l   

// If a is chosen from row[0] and e is chosen from row[1] then,  
float subCost1 = cost(a,e, path_to_a);  

// If i is chosen from row[2] ,  
float subCost2 = cost(e,i,path_to_e);  

// Cost for selecting aei combination is  
float cost1 = subCost1 + subCost2;  

//similarly other three costs need to be calculated by selecting other remaining elements  
//The elements should not intersect with each other eg. combinations aei and bej cannot exist on the same set.  

//Suppose the other combinations chosen are bfj with cost cost2, cgk with cost cost3   and dhl with cost cost4  
float totalCost = cost1 + cost2 + cost3 + cost4;   

//This is the cost got from one combination. All the other possible combinations should be enumerated to get the minimum cost combination. 

【问题讨论】:

如果你只是在寻找 C++ 的答案,正如标题所暗示的那样,为什么还要用 C 来标记它? 这是作业题吗?如果是这样,它应该被标记为这样 如果它们被存储为向量,C 标签可能不合适 这似乎是一个基于成本的图路由问题。如果您编辑以显示您已经拥有的一些代码,我们可能会为您提供一些开始。 @awoodland 好吧,如果有 c 算法,它也可以应用于 c++。所以我最初标记了 C。 【参考方案1】:

发布更多实用代码

见github:https://gist.github.com/1233012#file_new.cpp

这基本上是一种基于更简单的方法生成所有可能排列的更好方法(因此我之前没有真正的理由发布它:就目前而言,它不会除了python代码之外的任何东西)。

我还是决定分享它,因为您也许可以从中获得一些利润,作为最终解决方案的基础。

专业版:

快得多 更智能的算法(利用 STL 和数学 :)) 指令优化 存储优化 通用问题模型 模型和算法思想可以作为正确算法的基础 为良好的 OpenMP 并行化(n-方式,用于 n 行)设计的基础(但未充实)

对比:

以牺牲灵活性为代价,代码效率更高:使用更循序渐进的 Python 方法,调整代码以构建有关约束和成本启发式的逻辑会容易得多

总而言之,我觉得我的 C++ 代码可能是一个大赢家 IFF 事实证明,考虑到成本函数,模拟退火是合适的;代码中采用的方法会给

一种高效的存储模型 一种生成随机/密切相关的新网格配置的高效方法 便捷的显示功能

强制(abritrary...)基准数据点(与 python 版本的比较:)

  a  b  c  d e
  f  g  h  i j
  k  l  m  n o
  p  q  r  s t

Result: 207360000

real  0m13.016s
user  0m13.000s
sys   0m0.010s

这是我们到现在为止的:

从描述中我收集到建议您有一个基本图表,例如

必须构造一条访问网格中所有节点的路径 (Hamiltonian cycle)。

额外的限制是后续节点必须取自下一个rank(a-d、e-h、i-l 是三个rank em>;一旦访问了最后一个 rank 中的节点,路径必须与第一个 rank

中的任何未访问节点继续

边缘是加权的,因为它们具有相关的成本。然而,权重函数对于图算法来说并不是传统的,因为成本取决于完整路径,而不仅仅是每条边的端点。

鉴于此,我相信我们处于“全覆盖”问题的领域(需要 A* 算法,最著名的是 Knuths Dancing Links 论文)。

具体来说如果没有进一步的信息(路径的等价性,成本函数的特定属性),获得满足约束的“最便宜”汉密尔顿路径的最知名算法将是

生成所有可能的这样的路径 计算每个的实际成本函数 选择成本最低的路径

这就是为什么我开始编写一个非常愚蠢的蛮力生成器,它可以在 NxM 的通用网格中生成所有可能的唯一路径。

宇宙的尽头

3×4 样本网格的输出为 4!3 = 13824 条可能的路径...外推到 6×48 列,得到 6!48 = 1.4×10137 种可能性。很明显 如果不进一步优化,问题就难以解决(NP Hard 之类的——我不记得很微妙的定义)。

运行时的爆炸声震耳欲聋:

3×4(实测)耗时约0.175s 4×5(实测)大约需要 6 分 5 秒(在快速机器上在 PyPy 1.6 下无输出运行) 5×6 大约需要 10 年零 9 个月以上...

在 48x6 时,我们会看到...什么... 8.3x10107(仔细阅读)

现场观看:http://ideone.com/YsVRE

无论如何,这里是 python 代码(所有预设为 2×3 网格)

#!/usr/bin/python
ROWS = 2
COLS = 3

## different cell representations
def cell(r,c): 
    ## exercise for the reader: _gues_ which of the following is the fastest
    ## ...
    ## then profile it :)
    index = COLS*(r) + c
    # return [ r,c ]
    # return ( r,c )
    # return index
    # return "(%i,%i)" % (r,c)

    def baseN(num,b,numerals="abcdefghijklmnopqrstuvwxyz"):
        return ((num == 0) and numerals[0]) or (baseN(num // b, b, numerals).lstrip(numerals[0]) + numerals[num % b])

    return baseN(index, 26)

ORIGIN = cell(0,0)

def debug(t): pass; #print t
def dump(grid): print("\n".join(map(str, grid)))

def print_path(path):
    ## Note: to 'normalize' to start at (1,1) node:
    # while ORIGIN != path[0]: path = path[1:] + path[:1] 
    print " -> ".join(map(str, path))

def bruteforce_hamiltonians(grid, whenfound):
    def inner(grid, whenfound, partial):

        cols = len(grid[-1]) # number of columns remaining in last rank
        if cols<1:
            # assert 1 == len(set([ len(r) for r in grid ])) # for debug only
            whenfound(partial)                             # disable when benchmarking
            pass
        else:
            #debug(" ------ cols: %i ------- " % cols)

            for i,rank in enumerate(grid):
                if len(rank)<cols: continue
                #debug("debug: %i, %s (partial: %s%s)" % (i,rank, "... " if len(partial)>3 else "", partial[-3:]))
                for ci,cell in enumerate(rank):
                    partial.append(cell)
                    grid[i] = rank[:ci]+rank[ci+1:] # modify grid in-place, keeps rank

                    inner(grid, whenfound, partial)

                    grid[i] = rank # restore in-place
                    partial.pop()
                break
        pass
    # start of recursion
    inner(grid, whenfound, [])

grid = [ [ cell(c,r) for r in range(COLS) ] for c in range(ROWS) ]

dump(grid)

bruteforce_hamiltonians(grid, print_path)

【讨论】:

漂亮的图表。你是怎么做的? @Mooing Duck:添加源here -- 我使用 vimdot 来起草它们 @sehe:是的,你画的图对应我的问题。但是我还没有创建这个图形结构,我只有向量的向量和计算成本的函数,就像我在代码中写的那样。目标是找到到所有叶子 i、j、k、l 的最小成本路径集 a、b、c、d,并且这些路径不应相交。在我正在处理的问题中,最多有 48 列和 6 行。 @coolcav:两个问题:(a) 是重要字母还是无意义的占位符(您可以使用矩阵坐标代替)? (b) by intersect,你的意思是连接边交叉还是你的意思是路径重叠,因为一个节点被多次访问? -- 哦,PS:(c) 还有前面贴的关于成本函数的问题 @sehe:在我的实现中,字母代表映射到相应向量的整数。与这些数字相关的对应向量用于计算成本。例如,“a”可以映射到 [5,2,1,4],“e”可以映射到 [4,3,2,1]。这两个向量用于计算连接 a 和 b 的边的成本。我希望这能回答 (a) 和 (c)。【参考方案2】:

首先,有一个非常有用的观察结果。

我认为 4!^3 结果并没有捕捉到 aei, bfj, cgk, dhl 和(例如) bfj, aei, cgk, dhl 具有 same 的事实成本。

这意味着我们只需要考虑以下形式的序列

 a??, b??, c??, d?? 

这种等式将不同案例的数量减少了 4!

另一方面,@sehe 有 3x4 得到 4!^3(我同意),所以同样 6x48 需要 48!^6。这些“只有”48!^5 是不同的。现在是 2.95 × 10^305。

使用 3x4 示例,这是一个给出某种答案的算法的开始。

Enumerate all the triplets and their costs. 
Pick the lowest cost triplet.
Remove all remaining triplets containing a letter from that triplet.
Now find the lowest cost triplet remaining.
And so on.

请注意,这不是完全详尽的搜索。

从这里我也看到,这仍然是很多计算。第一次通过仍然需要计算 48^6 (12,230, 590, 464) 成本。我想这是可以做到的,但需要付出很多努力。相比之下,随后的通行证会很便宜。

【讨论】:

你和我的想法完全一样。但是,我打算询问 OP,成本函数是否实际上仅在单独的“vericals”上运行。 +1 在我之前发布这个想法:) @sehe。我正在根据 OP 评论“两种组合的路径 cgk 的成本将相同” +2 这个答案会让我投票赞成赏金:我还没有找到时间做任何腿部工作(如果我实际上相信这些假设,我可能会更有动力持有;我认为未来的程序员在开始编码之前必须与问题所有者更密切地合作)【参考方案3】:

编辑:添加完整的解决方案

正如其他答案已经指出的那样,您的问题太难以面对蛮力。此类问题的出发点始终是Simulated annealing。我创建了一个实现该算法的小应用程序。

查看问题的另一种方法是最小化复杂函数。此外,您对可能的解决方案有额外的限制。我从一个随机有效配置(满足您的约束)开始,然后我修改了该随机解决方案,每次更改一个元素。我强制应用程序执行有效的转换。代码很清晰。

我已经创建了一个模板函数,所以你只需要提供必要的函数对象和结构。

#include <iostream>
#include <cmath>
#include <ctime>
#include <vector>
#include <algorithm>
#include <functional>

//row [0]:  c00  c01  c02  c03   
//row [1]:  c10  c11  c12  c13  
//row [2]:  c20  c21  c22  c23 


typedef std::pair<int,int> RowColIndex;
// the deeper vector has size 3 (aei for example)
// the outer vector has size 4
typedef std::vector<std::vector<RowColIndex> > Matrix;

size_t getRandomNumber(size_t up)

    return rand() % up;


struct Configuration

    Configuration(const Matrix& matrix) : matrix_(matrix)
    Matrix matrix_;
;

std::ostream& operator<<(std::ostream& os,const Configuration& toPrint)

    for (size_t row = 0; row < toPrint.matrix_.at(0).size(); row++)
    
        for (size_t col = 0; col < toPrint.matrix_.size(); col++)
        
            os << toPrint.matrix_.at(col).at(row).first  << "," 
               << toPrint.matrix_.at(col).at(row).second << '\t';
        
        os << '\n';
       
    return os;


struct Energy 
 
    double operator()(const Configuration& conf)
    
        double result = 0;
        for (size_t col = 0; col < conf.matrix_.size(); col++)
        
            for (size_t row =0; row < conf.matrix_.at(col).size(); row++)
            
                result += pow(static_cast<double>(row) - static_cast<double>(conf.matrix_.at(col).at(row).first),2) +
                          pow(static_cast<double>(col) - static_cast<double>(conf.matrix_.at(col).at(row).second),2);
            
        
        return result;
    
;

size_t calculateNewColumn(std::vector<int>& isAlreadyUse)

    size_t random;
    do
    
        random = getRandomNumber(isAlreadyUse.size());
    
    while (isAlreadyUse.at(random) != 0);

    isAlreadyUse.at(random) = 1;
    return random;


Configuration createConfiguration(size_t numberOfRow,size_t numberOfColumn)

    //create suitable matrix
    Matrix matrix;
    //add empty column vector
    for (size_t col = 0; col < numberOfColumn; col++)
        matrix.push_back(std::vector<RowColIndex>());

    //loop over all the element
    for (size_t row = 0; row < numberOfRow; row++)
    
        std::vector<int> isAlreadyUse(numberOfColumn);
        for (size_t col = 0; col < numberOfColumn; col++)
        
            size_t newCol = calculateNewColumn(isAlreadyUse);
            matrix.at(newCol).push_back(std::make_pair(row,col));
        
       

    return Configuration(matrix);



struct CreateNewConfiguration

    Configuration operator()(const Configuration& conf)
    
        Configuration result(conf);

        size_t fromRow = getRandomNumber(result.matrix_.at(0).size());

        size_t fromCol = getRandomNumber(result.matrix_.size());
        size_t toCol = getRandomNumber(result.matrix_.size());

        result.matrix_.at(fromCol).at(fromRow) = conf.matrix_.at(toCol).at(fromRow);
        result.matrix_.at(toCol).at(fromRow) = conf.matrix_.at(fromCol).at(fromRow);

        return result;
    
;

template<typename Conf,typename CalcEnergy,typename CreateRandomConf>
Conf Annealing(const Conf& start,CalcEnergy energy,CreateRandomConf createNewConfiguration,
               int maxIter = 100000,double minimumEnergy = 1.0e-005)

    Conf first(start);
    int iter = 0;
    while (iter < maxIter && energy(first) > minimumEnergy )
    
        Configuration newConf(createNewConfiguration(first));
        if( energy(first) > energy(newConf))
        
            first = newConf;
        
        iter++;
    
    return first;


int main(int argc,char* argv[])


    size_t nRows = 25;
    size_t nCols = 25;
    std::vector<Configuration> res;
    for (int i =0; i < 10; i++)
    
        std::cout << "Configuration #" << i << std::endl;
        Configuration c = createConfiguration(nRows,nCols);
        res.push_back(Annealing(c,Energy(),CreateNewConfiguration()));
    

    std::vector<Configuration>::iterator it = res.begin();


    std::vector<Configuration>::iterator lowest = it;
    while (++it != res.end())
    
        if (Energy()(*it) < Energy()(*lowest))
            lowest = it;
    

    std::cout << Energy()(*lowest) << std::endl;

    std::cout << std::endl;

    std::cout << *lowest << std::endl;


    std::cin.get();
    return 0;

当然,您不能保证解决方案是最好的(这是一种启发式方法)。不过,这是一个很好的起点。

您还没有提供完整的功能成本,所以我实现了自己的成本,让我可以简单地检查最终结果。您只需要提供功能成本即可完成。

你可能会让程序更高效,还有很大的改进空间,但逻辑是存在的,你可以很容易地实现你的功能。

复杂性

算法的复杂度为 E*I*C 其中 I = 迭代次数 C = 随机配置的数量(避免局部最小值) E = 计算能量函数(或函数成本)

在这种情况下,E 实际上是 N*M,其中 N 和 M 是初始矩阵的维数。

如果您对模拟退火结果不满意,可以尝试genetic algorithms。

【讨论】:

+1,绝对使用模拟退火或遗传算法解决这个问题。 到目前为止,我认为示例并不简单,也不完整。成本函数至今未指定。我不能假设它适合 SA。 如果是这样,这是一个很好的选择(实施可能会更有效)。但是现在请参阅WP:问题在于确实没有任何逻辑neighbour 生成器(能量函数在随机选择的转换处可能是离散的甚至是奇异的。同样,从问题域应用启发式方法可以完全连接这些差距,如果他们知道的话 我尽我所能利用可用的信息量,我知道实施效率不高,但我认为清晰总比快速好。代码的目的是展示一种算法。【参考方案4】:

你可以递归解决问题。

该方法的输入是要计算的第一个向量的索引,该向量在函数外部共享。

对于剩余两行的情况,可以使用回溯计算解决方案。在这种情况下,您只需要找到更便宜的配对即可。

对于多于两行的情况,应该调用下一个索引的方法,获取部分结果,然后再次使用回溯计算最小值。

当流回到第一个向量时,可以将结果合并为最终结果。

【讨论】:

【参考方案5】:

值得注意的是,对于一些有趣的路径成本选择,有一个多时间算法,例如,如果路径成本是边缘成本的总和,那么可以通过运行找到最优解,对于所有 i,第 i 行和第 i + 1 行的匈牙利算法。

【讨论】:

以上是关于在 C++ 中通过网格/矩阵找到成本优化路径的主要内容,如果未能解决你的问题,请参考以下文章

移动度数为 8 时矩阵中的最小成本路径

使用向右或向下跳跃或单位步长在网格中的最小成本路径

尝试在 C++ 中通过引用乘以矩阵时出现访问冲突异常

如何优化软件测试成本?——转译

如何在 Floyd-Warshall 算法中找到最短路径和最短成本

如何找到成本最高的路径