使用pyCUDA实现批量化的Homotopy Continuation求解器

使用pyCUDA实现批量化的Homotopy Continuation求解器

GPU-HC

Homotopy Continuation可以用来求解多项式方程，原理可以看之前的一篇论文阅读--在RANSAC下基于HC和ML的位姿估计最小配置问题。而对于PnP问题而言，由于可以把位姿估计问题转化成求解一系列方程组的问题，且这些方程组天然存在多个伪解，所以使用HC去求解也同样具有优势。

由于上述过程相当于每个HC跟踪一个解，在Minimal Solutions to Generalized Three-View Relative Pose Problem和GPU-Based Homotopy Continuation for Minimal Problems in Computer Vision这两篇文章中，就使用CUDA来编写了HC求解器，使用Monodromy作为start system，每个CUDA core跟踪一个HC曲线，利用Monodromy的性质，先对具体的问题找出所有的参数，然后每个GPU-HC就可以利用这些参数来针对一个具体的问题求解出所有的解。关于Monodromy、Polyhedral这些HC start system可以简单参考这里，上面两篇论文的Monodromy参数计算也是基于这个库，这个库是julia编写的，可以直接用jupyter运行。

python中的HC

但是上述代码是纯CUDA代码，而且还需要使用julia编写的程序计算monodromy参数，如果想要在python，尤其是和pytorch中结合使用会稍微有些麻烦。而且HC作为一个优化步骤，涉及到的循环很多，如果用pytorch的autogradient来写的话效率其实很差，并且如果把HC作为网络层后面的优化器的话，那么纯python实现必须把整个batch的数据展开成一维来逐个处理，计算量爆炸。

所以基于此，可以使用pyCUDA来实现一个batch_HC，这样就能跟torch的网络层无缝衔接。

简单测试一下性能，每个HC的超参数为：t_size=0.05 newton_iter=400(此处代码暂不公开)

数据量	python-HC	pyCUDA-HC
1	5.05s	4.72+0.032s
5	24.99s	4.51+0.027s
10	50.61s	5.15+0.033s
100	未测	5.04+0.031s

由于pyCUDA其实就是先用CUDA写一遍，然后再把CUDA代码转化成python接口，所以本质上还是在用C++写CUDA，因此核函数在pyCUDA中也同样需要编译，所以pyCUDA分为编译时间(前)和计算时间(后)，在我个人电脑(rtx2060，CUDA10.2)编译耗时大概为5s左右，而计算时间则是0.03s左右，且无论针对多么大的batch，运行时间基本都一致；而使用torch实现的HC，基本上稳定计算一个数据为5s左右，而通过torch写的代码很难向量化，所以导致运行效率低下。