CN108053065A - 一种基于gpu绘制的半离散最优传输方法及系统 - Google Patents
一种基于gpu绘制的半离散最优传输方法及系统 Download PDFInfo
- Publication number
- CN108053065A CN108053065A CN201711311826.XA CN201711311826A CN108053065A CN 108053065 A CN108053065 A CN 108053065A CN 201711311826 A CN201711311826 A CN 201711311826A CN 108053065 A CN108053065 A CN 108053065A
- Authority
- CN
- China
- Prior art keywords
- discrete
- gpu
- optimal transmission
- source region
- measure
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q10/00—Administration; Management
- G06Q10/04—Forecasting or optimisation specially adapted for administrative or management purposes, e.g. linear programming or "cutting stock problem"
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Business, Economics & Management (AREA)
- Mathematical Physics (AREA)
- Theoretical Computer Science (AREA)
- Pure & Applied Mathematics (AREA)
- Human Resources & Organizations (AREA)
- Data Mining & Analysis (AREA)
- Economics (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- Strategic Management (AREA)
- Computational Mathematics (AREA)
- Computing Systems (AREA)
- Development Economics (AREA)
- General Business, Economics & Management (AREA)
- Entrepreneurship & Innovation (AREA)
- Tourism & Hospitality (AREA)
- Quality & Reliability (AREA)
- Operations Research (AREA)
- Game Theory and Decision Science (AREA)
- Algebra (AREA)
- Marketing (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Image Generation (AREA)
- Medicines That Contain Protein Lipid Enzymes And Other Medicines (AREA)
Abstract
本发明公开了一种基于GPU绘制的半离散最优传输方法及系统，首先输入连续的源区域，并设置绘制视口大小；输入离散母点，给每个离散母点赋予一个狄拉克测度和权重；生成随机颜色表；每次迭代过程中，对每个离散母点，构造空间四边形；利用GPU对空间四边形进行渲染；求解空间四边形的上包络，并将其向源区域正投影得到源区域的胞腔分解；计算每个加权维诺图胞腔的测度；计算传输映射的能量函数及其梯度；更新离散母点的权重；循环直到能量函数达到极小值，得到最优传输的最终解。本发明原理简单易懂，实现简单高效，避免了复杂的胞腔分解、积分计算，克服了CPU计算效率受限的问题，大大提高了最优传输问题的求解效率，适用于求解大规模最优传输问题。
Description
技术领域
本发明属于离散最优传输技术领域，具体应用广泛，包括三维模型处理，如三维网格参数化、曲面注册等；包括计算机图形图像处理，如颜色传输等；包括深度学习，如生成对抗网络WGAN等，涉及一种基于GPU绘制的半离散最优传输方法及系统。
背景技术
1781年，法国数学家蒙日最早提出了最优传输问题：假设X,Y是完备可分的欧式区域，其上分别定义概率测度μ,ν；传输代价函数c(x,y):X×Y→[0,+∞]为光滑函数，表示质量点x∈X传输到点y∈Y的传输代价；从欧式区域X到欧式区域Y的映射T:X→Y将概率测度μ前推到欧式区域Y上的概率测度T#μ。如果对于任意的波莱尔集都满足T#μ(B)＝μ(T-1(B))，则称映射T为保测度映射，所有保测度映射中传输代价最小的映射为最优传输映射：inf{M(T):＝∫Xc(x,T(x))dμ(x):T#μ(B)＝μ}。蒙日提出的传输代价为L1距离，即c(x,y)＝|x-y|，数学分析困难。1948年，俄罗斯数学家康塔洛维奇运用线性规划方法，将传输映射推广为传输方案，求解了最优传输问题，并证明了最优传输问题解的存在性和唯一性。最优传输问题带来巨大经济价值，康塔洛维奇因此获得了1975年的诺贝尔经济学奖。康塔洛维奇的线性规划思想是：将欧式区域X,Y离散化，表示成一簇离散点X＝{xi}和Y＝{yj}，概率测度μ,ν离散化成狄拉克测度和最优传输方案γ离散化成狄拉克测度这样，最优传输方案问题被转换成经典的线性规划问题：s.t.∑iγij＝νj,∑jγij＝μi。1981年，法国数学家Brenier证明：当传输代价c(x,y)＝|x-y|2，且欧式空间X为凸区域时，某个凸函数f的梯度即为最优传输映射，建立了最优传输映射理论和凸几何之间的内在联系。顾险峰等结合计算几何与凸优化的方法，对Brenier理论进行离散化，建立了欧式空间中凸区域上离散最优传输的凸几何理论，构造出最优传输问题的离散求解算法。
上述求解最优传输问题的方法各具优缺点，康塔洛维奇提出的线性规划方法原理简单，但实现困难且复杂度高，空间复杂度为N2(N为离散点的个数)，计算效率低，无法适用于大规模的最优传输问题。例如，三维模型往往具有上万个顶点，基于线性规划方法的最优传输问题求解方案在三维模型的处理中，效率低下。顾险峰提出的计算几何凸方法将求解问题的规模从N2降到N，极大地降低了空间复杂度，提高了计算效率。但目前的最优传输求解方法几乎全是基于CPU计算的，计算效率受限。
对最优传输问题的研究在全世界范围内得到重视，成果丰硕，极大地推动了数学、经济学、医学、人工智能、计算机科学、医学等学科的发展。这些学科的高速发展反过来推动最优传输问题往更高求解效率的方向发展，但受CPU计算效率限制，目前最优传输求解方法，效率有待进一步提高，以适用于规模更大、效率要求更高的最优传输问题中。
发明内容
为了解决上述技术问题，本发明提供了一种基于GPU绘制的半离散最优传输方法及系统，能够极大地提高最优传输问题的求解效率。
本发明的方法所采用的技术方案是：一种基于GPU绘制的半离散最优传输方法，其特征在于，包括以下步骤：
步骤1：输入连续的源区域，并设置绘制视口大小；输入离散母点，给每个离散母点赋予一个狄拉克测度和权重；
步骤2：生成随机颜色表；
步骤3：对每个离散母点，构建空间四边形；
步骤4：利用GPU绘制空间四边形；
步骤5：构造空间四边形的上包络，并将其向源区域正投影得到源区域的胞腔分解；计算每个加权维诺图胞腔的测度；
步骤6：计算传输映射的能量函数及其梯度；
步骤7：更新离散母点的权重；
步骤8：循环步骤3-步骤7，直到能量函数达到极小值，得到最优传输的最终解。
本发明的系统所采用的技术方案是：一种基于GPU绘制的半离散最优传输系统，其特征在于：包括以下模块；
模块1：用于输入连续的源区域，并设置绘制视口大小；输入离散母点，给每个离散母点赋予一个狄拉克测度和权重；
模块2：用于生成随机颜色表；
模块3：用于对每个离散母点，构建空间四边形；
模块4：用于利用GPU绘制空间四边形；
模块5：用于构造空间四边形的上包络，并将其向源区域正投影得到源区域的胞腔分解；计算每个加权维诺图胞腔的测度；
模块6：用于计算传输映射的能量函数及其梯度；
模块7：用于更新离散母点的权重；
模块8：用于循环步骤3-步骤7，直到能量函数达到极小值，得到最优传输的最终解。
本发明的有益效果在于：本发明提供了一种基于GPU绘制的半离散最优传输方法，将目标区域离散化，表示成一簇离散母点。对每个离散母点，构造定义在源区域上的空间四边形。用户只需定义空间四边形的顶点坐标，GPU自动完成绘制。GPU自动根据像素点DepthBuffer的最大值求解空间四边形的上包络，简单高效，无需用户干预。空间四边形上包络向绘制视口正投影获得源区域的胞腔分解。胞腔的测度通过Compute Shader快速统计像素点个数的方法求得，避免了复杂的多边形分解及积分计算。本发明提供的基于GPU绘制的半离散最优传输方法，原理简单易懂，实现简单高效，避免了复杂的胞腔分解、积分计算，克服了CPU计算效率受限的问题，大大提高了最优传输问题的求解效率，适用于求解大规模最优传输问题。
附图说明
图1为本发明实施例的半离散最优传输问题示意图；
图2为本发明实施例的半离散最优传输问题的求解原理图；
图3为本发明实施例的基于GPU绘制的半离散最优传输方法流程图；
图4为本发明实施例的GPU绘制空间四边形流程图；
图5为本发明实施例的离散母点狄拉克测度相等的最优传输结果示例图；
图6为本发明实施例的离散母点狄拉克测度不相等的最优传输结果示例图；
图7为本发明实施例的基于半离散最优传输的网格保面积参数化示例图。
具体实施方式
为了便于本领域普通技术人员理解和实施本发明，下面结合附图及实施例对本发明作进一步的详细描述，应当理解，此处所描述的实施示例仅用于说明和解释本发明，并不用于限定本发明。
请见图1，本发明的半离散最优传输问题是：源区域为连续区域，概率测度分布均匀，测度密度函数处处为μ；目标区域被离散化，表示成一簇离散母点，每个离散母点yj上被赋予一个狄拉克测度νj。求解将源区域映射到目标区域离散母点的最优传输映射f，使得映射f满足以下条件：映射前后测度相等，即且传输代价最小，即映射的能量取得极小值。
请见图2，本发明的半离散最优传输问题的求解原理为：给定截距向量，构造定义在源区域上的超平面，超平面的上包络为分段线性凸函数的图；上包络在源区域的正向投影诱导源区域的胞腔分解，线性凸函数的梯度将每个加权维诺图胞腔映射到相应的离散母点。
为了解决本发明提出的半离散最优传输问题，请见图3，本发明提供了基于GPU绘制的半离散最优传输方法，包括以下步骤：
步骤1：输入连续的源区域，并设置绘制视口大小；输入离散母点，给每个离散母点赋予一个狄拉克测度和权重；
在本实施方式中，输入的源区域X为连续的二维欧式区域，概率测度分布均匀，测度的密度函数为μ(x)。目标区域Y被离散化，表示为一簇离散母点Y＝{y1,y2,…,yM}；离散母点yj的坐标为yj＝(Aj,Bj),Aj,Bj的值为(-1,1)间的随机数，母点yj上赋予的狄拉克测度为νj，权重为hj。本实施方式中，权重hj初始化母点到原点的距离平方的一半，即其中<,>表示点积。源区域的测度与目标区域的狄拉克测度相等，即满足以下条件：
其中，M为离散母点的总数。
将绘制视口设置为W×H的平面四边形(为母点构造空间四边形，绘制的视口是空间母点的四边形在源区域的正投影，故为平面四边形)，W表示宽度，H表示高度，绘制视口的四个顶点坐标分别为：
其中，
步骤2：生成随机颜色表；
GPU仅支持绘制三角面片，因此将一个空间四边形分解为两个三角形。根据母点个数M，生成大小为2M的随机颜色表Colors，每种颜色由r,g,b,a四个通道组成。通道r,g,b的值为(0,1)之间的随机数，a通道存储颜色编号。属于同一个母点yj所在四边形的两个三角面片对应相同的颜色，Colors[j+1]＝Colors[j]。
步骤3：对每个离散母点，构造空间四边形；
对离散母点yj，构造空间四边形Qj，其顶点存储绘制视口的顶点坐标信息、母点yj的坐标信息和权重，以及颜色编号，具体如下：
其中
yj＝(Aj,Bj)、hj分别表示母点yj的顶点坐标和权重。
步骤4：利用GPU对空间四边形进行渲染；
请见图4，利用GPU绘制空间四边形的具体步骤如下：
步骤4.1：将空间四边形Qj划分为两个三角面片，Tj1＝(Vj0,Vj1,Vj2)和Tj2＝(Vj2,Vj1,Vj3)，将空间四边形的顶点信息Vertices，三角面片的拓扑结构T传入Vertex Shader中；
步骤4.2：Vertex Shader根据三角面片顶点信息中的颜色编号在随机颜色表Colors中找到对应的颜色，对顶点进行着色；
步骤4.3：GPU对三角面片进行光栅化处理得到Fragment；
步骤4.4：Pixel/Fragment Shader对Fragment进行着色渲染；
步骤4.5：绘制图元，显示在绘制视口界面上。
步骤5：求解空间四边形的上包络，并将其向源区域正投影得到源区域的胞腔分解；计算每个加权维诺图胞腔的测度；
定义源区域X上的凸线性函数：
其中，M表示离散母点的个数。线性函数fh(x)存储在空间四边顶点的Z轴上。fh(x)的图即为步骤4所得的空间四边形的上包络ε(h)。GPU根据像素点Depth Buffer的最大值构造空间四边形的上包络。将上包络ε(h)向源区域正投影，得到源区域的胞腔分解υ(h)，即加权维诺图胞腔为Wj(h)＝{x∈X|μh(x)＝<x,yj>+hj}∩X，其测度计算公式如下：
对于概率测度分布均匀的源区域，通过GPU的Compute Shader计算每片Fragment的像素点个数，根据a通道的颜色编号统计相同颜色编号的像素点总个数，以获得维诺图胞腔Wj(h)的测度wj(h)。
步骤6：计算传输映射的能量函数及其梯度；
迭代过程中，按照以下公式计算半离散最优传输的能量函数，
其中，fh(x)为凸线性函数，μ(x)为源区域的密度函数；h表示截距向量，hj为离散母点yj的权重，vj为离散母点yj上赋予的狄拉克测度，M为离散母点的个数。
能量函数梯度的计算公式如下：
步骤7：更新离散母点的权重；
离散母点权重的更新规则如下：
其中，stepLength为下降步长，步长越大，能量函数收敛越快，但结果精确度越低；步长越小，结果精确度越高，但收敛速度越慢。可采用拟牛顿下降法等方法来更新截距向量h；为了获得较快的收敛速度和较精确的结果，本实施例中，下降步长stepLength设置为0.001；表示能量函数在胞腔Wj(h)的梯度。
步骤8：循环步骤3-步骤7，直到能量函数达到极小值，得到最优传输的最终解；
母点的权重更新后，空间四边形也更新，由空间四边形的上包络诱导的胞腔分解也得到更新。因此，每个胞腔的测度发生改变，当每个胞腔的测度与母点上赋予的狄拉克测度相等时，能量函数达到极小值，算法终止，得到最终的最优传输结果。
能量函数达到极小值的判断条件是，其梯度函数收敛，即梯度函数的最大值小于预设的收敛误差，
max(wj(h)-νj)＜threshold
其中，wj(h)为加权维诺图胞腔Wj(h)的测度，νj为离散母点上赋予的狄拉克测度。在本实施例中，收敛阈值设为10-6。
在本实施例中，通过将目标区域离散化表示成离散母点，对每个离散母点，构造定义在源区域的空间四边形。将复杂的上包络求解问题转换成GPU求解像素点Depth Buffer最大值问题，将复杂的多边形分解，计算积分转成简单的统计编号相同的像素点个数问题，极大地提高了计算效率。本发明提供的基于GPU绘制的半离散最优传输方法，原理简单易懂，实现简单高效，避免了复杂的胞腔分解、积分计算，克服了CPU计算效率受限的问题，大大提高了最优传输问题的求解效率，适用于求解大规模最优传输问题。图5为离散母点狄拉克测度相等的最优传输结果示例，其中，源区域为大小800×800的平面四边形，概率测度均匀分布，总测度为∫Xdμ；绘制视口的大小与源区域相同。目标区域离散化，表示为一簇离散母点，其上赋予的狄拉克测度为M＝1000为离散母点的总个数，母点权重初始化为母点到原点距离平方的一半。图(a)为原始胞腔分解；图(b)为中间结果，此时，存在某个胞腔Wj(h)，其面积wj(h)与该母点上的狄拉克测度νj不等；图(c)为最终结果，此时，所有胞腔的面积与母点上给定的狄拉克测度相等。当离散母点的狄拉克测度相等时，所有胞腔的面积都相等。图6为离散母点狄拉克测度不相等的最优传输结果示例，当离散母点的狄拉克测度不相等时，源区域胞腔的面积也不相等。图7为基于本发明提出的基于GUP绘制的半离散最优传输方法的网格保面积参数化。图(a)为原始三维网格模型；图(b)是原始三维网格模型通过曲面RICCI曲率流保角参数化到单位圆盘，该参数化诱导较大的面积畸变；图(c)是利用本发明提出的基于GUP绘制的半离散最优传输方法对图(b)参数化的面积畸变进行修正后的保面积参数化。其中，图(b)保角参数化为源区域，离散母点为图(b)三维网格的顶点，狄拉克测度为图(a)每个顶点一周三角形面积之和的三分之一。绘制视口为800×800的平面四边形。
本发明还提供了一种基于GPU绘制的半离散最优传输系统，包括以下模块；
模块1：用于输入连续的源区域，并设置绘制视口大小；输入离散母点，给每个离散母点赋予一个狄拉克测度和权重；
模块2：用于生成随机颜色表；
模块3：用于对每个离散母点，构建空间四边形；
模块4：用于利用GPU绘制空间四边形；
模块5：用于构造空间四边形的上包络，并将其向源区域正投影得到源区域的胞腔分解；计算每个加权维诺图胞腔的测度；
模块6：用于计算传输映射的能量函数及其梯度；
模块7：用于更新离散母点的权重；
模块8：用于循环步骤3-步骤7，直到能量函数达到极小值，得到最优传输的最终解。
以上优选实施例虽然只是二维的，但是本发明提出的基于GPU绘制的半离散最优传输方法同样适用于求解三维半最优传输问题，技巧是将三维区域按Z轴进行分片处理。除了优选实施例中的网格保面积参数化，本发明提出的基于GPU绘制的半离散最优传输方法同样适用于其他最优传输问题，如曲面注册、颜色转移以及生成对抗网络WGAN等技术领域。
最后说明的是，以上优选实施例仅用以说明本发明的技术方案而非限制，尽管通过上述优选实施例已经对本发明进行了详细的描述，但本领域技术人员应当理解，可以在形式上和细节上对其作出各种各样的改变，而不偏离本发明权利要求书所限定的范围。
Claims (10)
1.一种基于GPU绘制的半离散最优传输方法，其特征在于，包括以下步骤：
步骤1：输入连续的源区域，并设置绘制视口大小；输入离散母点，给每个离散母点赋予一个狄拉克测度和权重；
步骤2：生成随机颜色表；
步骤3：对每个离散母点，构建空间四边形；
步骤4：利用GPU绘制空间四边形；
步骤5：构造空间四边形的上包络，并将其向源区域正投影得到源区域的胞腔分解；计算每个加权维诺图胞腔的测度；
步骤6：计算传输映射的能量函数及其梯度；
步骤7：更新离散母点的权重；
步骤8：循环步骤3-步骤7，直到能量函数达到极小值，得到最优传输的最终解。
2.根据权利要求1所述的基于GPU绘制的半离散最优传输方法，其特征在于：步骤1中，源区域为连续的欧式区域X，概率测度分布均匀，处处为μ，绘制视口设置为大小W×H的平面四边形，W表示宽度，H表示高度，绘制视口的四个顶点坐标分别为：
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>p</mi>
<mn>0</mn>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>X</mi>
<mo>,</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>Y</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>p</mi>
<mn>1</mn>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>X</mi>
<mo>,</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>Y</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>p</mi>
<mn>2</mn>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>X</mi>
<mo>,</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>Y</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>p</mi>
<mn>3</mn>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>X</mi>
<mo>,</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>Y</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
其中，
目标区域Y离散化，表示成一簇离散母点Y＝(y1,y2,…,yM)，母点yj的坐标为(Aj,Bj)，其上赋予的狄拉克测度为νj，权重为hj；目标区域的总测度与源区域的测度相等，即满足以下条件：
<mrow>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>M</mi>
</munderover>
<msub>
<mi>v</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<msub>
<mo>&Integral;</mo>
<mi>X</mi>
</msub>
<mi>d</mi>
<mi>&mu;</mi>
</mrow>
其中，M表示离散母点的个数。
3.根据权利要求2所述的基于GPU绘制的半离散最优传输方法，其特征在于：步骤2中，颜色表Colors的大小为2M，其中M表示离散母点的个数；颜色表的颜色值由四个通道组成，即Colors[j]＝(r,g,b,a)，通道r,g,b的值为(0,1)之间的随机数，a通道存储颜色编号；属于同一个母点yj所在空间四边形的两个三角面片对应相同的颜色，Colors[j+1]＝Colors[j]。
4.根据权利要求2所述的基于GPU绘制的半离散最优传输方法，其特征在于：步骤3中，对每个离散母点yj，构造一个空间四边形Qj，其顶点存储绘制视口的顶点坐标信息、母点yj的信息、以及颜色编号，具体如下：
<mrow>
<msub>
<mi>Vertices</mi>
<msub>
<mi>Q</mi>
<mi>j</mi>
</msub>
</msub>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>V</mi>
<mrow>
<mi>j</mi>
<mn>0</mn>
</mrow>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>LB</mi>
<mi>j</mi>
</msub>
<mo>,</mo>
<msub>
<mi>CID</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>V</mi>
<mrow>
<mi>j</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mn>1</mn>
</msub>
<mo>,</mo>
<msub>
<mi>LT</mi>
<mi>j</mi>
</msub>
<mo>,</mo>
<msub>
<mi>CID</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>V</mi>
<mrow>
<mi>j</mi>
<mn>2</mn>
</mrow>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mn>2</mn>
</msub>
<mo>,</mo>
<msub>
<mi>RB</mi>
<mi>j</mi>
</msub>
<mo>,</mo>
<msub>
<mi>CID</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>V</mi>
<mrow>
<mi>j</mi>
<mn>3</mn>
</mrow>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<msub>
<mi>p</mi>
<mn>3</mn>
</msub>
<mo>,</mo>
<msub>
<mi>RT</mi>
<mi>j</mi>
</msub>
<msub>
<mi>CID</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
</mrow>
其中
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>LB</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mo><</mo>
<msub>
<mi>p</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>y</mi>
<mi>j</mi>
</msub>
<mo>></mo>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>X</mi>
<mo>&times;</mo>
<msub>
<mi>A</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>Y</mi>
<mo>&times;</mo>
<msub>
<mi>B</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>LT</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mo><</mo>
<msub>
<mi>p</mi>
<mn>1</mn>
</msub>
<mo>,</mo>
<msub>
<mi>y</mi>
<mi>j</mi>
</msub>
<mo>></mo>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>X</mi>
<mo>&times;</mo>
<msub>
<mi>A</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>Y</mi>
<mo>&times;</mo>
<msub>
<mi>B</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>RB</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mo><</mo>
<msub>
<mi>p</mi>
<mn>2</mn>
</msub>
<mo>,</mo>
<msub>
<mi>y</mi>
<mi>j</mi>
</msub>
<mo>></mo>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>X</mi>
<mo>&times;</mo>
<msub>
<mi>A</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mi>Y</mi>
<mo>&times;</mo>
<msub>
<mi>B</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>RB</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mo><</mo>
<msub>
<mi>p</mi>
<mn>3</mn>
</msub>
<mo>,</mo>
<msub>
<mi>y</mi>
<mi>j</mi>
</msub>
<mo>></mo>
<mo>+</mo>
<msub>
<mi>C</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>X</mi>
<mo>&times;</mo>
<msub>
<mi>A</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<mi>M</mi>
<mi>a</mi>
<mi>x</mi>
<mi>Y</mi>
<mo>&times;</mo>
<msub>
<mi>B</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
yj＝(Aj,Bj)、hj分别表示母点yj的顶点坐标和权重。
5.根据权利要求4所述的基于GPU绘制的半离散最优传输方法，其特征在于，步骤4的具体实现包括以下子步骤：
步骤4.1：将空间四边形Qj划分为两个三角面片，Tj1＝(Vj0,Vj1,Vj2)和Tj2＝(Vj2,Vj1,Vj3)，将空间四边形的顶点信息Vertices Qj、三角面片的拓扑结构T传入Vertex Shader；
步骤4.2：Vertex Shader根据顶点信息中的颜色编号在随机颜色表Colors中找到对应的颜色，对顶点进行着色；
步骤4.3：GPU对三角面片进行光栅化处理得到Fragment；
步骤4.4：Pixel/Fragment Shader对Fragment进行着色渲染；
步骤4.5：绘制图元，显示在绘制视口界面上。
6.根据权利要求4所述的基于GPU绘制的半离散最优传输方法，其特征在于：步骤5中，定义源区域X上的凸线性函数：
<mrow>
<msub>
<mi>f</mi>
<mi>h</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>:</mo>
<mo>=</mo>
<munder>
<mi>max</mi>
<mrow>
<mn>1</mn>
<mo>&le;</mo>
<mi>j</mi>
<mo>&le;</mo>
<mi>M</mi>
</mrow>
</munder>
<mo>{</mo>
<mo><</mo>
<mi>x</mi>
<mo>,</mo>
<msub>
<mi>y</mi>
<mi>j</mi>
</msub>
<mo>></mo>
<mo>+</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
<mo>}</mo>
</mrow>
其中，M表示离散母点的个数；fh(x)的图即为空间四边形的上包络ε(h)；将空间四边形的上包络ε(h)向源区域正投影，得到源区域的胞腔分解υ(h)，即
加权维诺图胞腔为Wj(h)＝{x∈X|fh(x)＝<x,yj>+hj}∩X，其测度计算公式如下：
<mrow>
<msub>
<mi>w</mi>
<mi>j</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mo>&Integral;</mo>
<mrow>
<msub>
<mi>W</mi>
<mi>j</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
</mrow>
</msub>
<mi>d</mi>
<mi>&mu;</mi>
</mrow>
GPU根据像素点Depth Buffer的最大值构造空间四边形的上包络；对于概率测度分布均匀的源区域，通过GPU的Compute Shader计算每片Fragment的像素点个数，根据a通道的颜色编号统计相同颜色编号的像素点总个数；由于颜色表的颜色编号与空间四边形的编号一致，颜色编号为j的像素点总数即为加权维诺图胞腔Wj(h)的测度wj(h)。
7.根据权利要求6所述的基于GPU绘制的半离散最优传输方法，其特征在于：步骤6中，传输映射的能量函数的计算公式如下：
<mrow>
<mi>E</mi>
<mrow>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
<mo>:</mo>
<mo>=</mo>
<msub>
<mo>&Integral;</mo>
<mi>X</mi>
</msub>
<msub>
<mi>f</mi>
<mi>h</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>&mu;</mi>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>M</mi>
</munderover>
<msub>
<mi>v</mi>
<mi>j</mi>
</msub>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
</mrow>
其中，fh(x)为凸线性函数，μ(x)为源区域的密度函数；h表示截距向量，hj为离散母点yj的权重，vj为离散母点yj上赋予的狄拉克测度，M为离散母点的个数；
能量函数梯度的计算公式如下：
<mrow>
<mo>&dtri;</mo>
<mi>E</mi>
<mrow>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mrow>
<mo>(</mo>
<msub>
<mi>w</mi>
<mn>1</mn>
</msub>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>v</mi>
<mn>1</mn>
</msub>
<mo>,</mo>
<msub>
<mi>w</mi>
<mn>2</mn>
</msub>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>v</mi>
<mn>2</mn>
</msub>
<mo>,</mo>
<mn>...</mn>
<mo>,</mo>
<msub>
<mi>w</mi>
<mi>M</mi>
</msub>
<mo>(</mo>
<mi>h</mi>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>v</mi>
<mi>M</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>.</mo>
</mrow>
8.根据权利要求7所述的基于GPU绘制的半离散最优传输方法，其特征在于：步骤7中，离散母点权重的更新规则如下：
<mrow>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
<mo>+</mo>
<mi>s</mi>
<mi>t</mi>
<mi>e</mi>
<mi>p</mi>
<mi>L</mi>
<mi>e</mi>
<mi>n</mi>
<mi>g</mi>
<mi>t</mi>
<mi>h</mi>
<mo>&times;</mo>
<mo>&dtri;</mo>
<mi>E</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>h</mi>
<mi>j</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
其中，stepLength为下降步长，步长越大，能量函数收敛越快，但结果精确度越低；步长越小，结果精确度越高，但收敛速度越慢；采用拟牛顿下降法法更新截距向量h；表示能量函数在胞腔Wj(h)的梯度。
9.根据权利要求7所述的基于GPU绘制的半离散最优传输方法，其特征在于：步骤8中，能量函数达到极小值的判断条件是，其梯度函数收敛，即梯度函数的最大值小于预设的收敛误差，
max(wj(h)-νj)＜threshold
其中，wj(h)为加权维诺图胞腔Wj(h)的测度，νj为离散母点上赋予的狄拉克测度。
10.一种基于GPU绘制的半离散最优传输系统，其特征在于：包括以下模块；
模块1：用于输入连续的源区域，并设置绘制视口大小；输入离散母点，给每个离散母点赋予一个狄拉克测度和权重；
模块2：用于生成随机颜色表；
模块3：用于对每个离散母点，构建空间四边形；
模块4：用于利用GPU绘制空间四边形；
模块5：用于构造空间四边形的上包络，并将其向源区域正投影得到源区域的胞腔分解；计算每个加权维诺图胞腔的测度；
模块6：用于计算传输映射的能量函数及其梯度；
模块7：用于更新离散母点的权重；
模块8：用于循环步骤3-步骤7，直到能量函数达到极小值，得到最优传输的最终解。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711311826.XA CN108053065B (zh) | 2017-12-11 | 2017-12-11 | 一种基于gpu绘制的半离散最优传输方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711311826.XA CN108053065B (zh) | 2017-12-11 | 2017-12-11 | 一种基于gpu绘制的半离散最优传输方法及系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN108053065A true CN108053065A (zh) | 2018-05-18 |
| CN108053065B CN108053065B (zh) | 2021-08-03 |
Family
ID=62124135
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201711311826.XA Active CN108053065B (zh) | 2017-12-11 | 2017-12-11 | 一种基于gpu绘制的半离散最优传输方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN108053065B (zh) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111008668A (zh) * | 2019-12-09 | 2020-04-14 | 上海熵熵微电子科技有限公司 | 一种数据奇异点的识别方法、装置、设备和介质 |
| CN111325337A (zh) * | 2019-08-19 | 2020-06-23 | 大连理工大学 | 一种基于自编码神经网络与最优传输的生成模型 |
| CN111950045A (zh) * | 2020-08-03 | 2020-11-17 | 大连理工大学 | 一种无光能损失的透镜设计方法 |
| CN113064272A (zh) * | 2021-03-04 | 2021-07-02 | 武汉大学 | 半离散最优传输下的光学自由曲面构造方法及系统 |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20060285114A1 (en) * | 2005-06-17 | 2006-12-21 | Canon Kabushiki Kaisha | Gas detection and photonic crystal devices design using predicted spectral responses |
| US20140192046A1 (en) * | 2013-01-07 | 2014-07-10 | Ecole Centrale Paris | Method and device for elastic registration between a two-dimensional digital image and a slice of a three-dimensional volume with overlapping content |
| US20150293884A1 (en) * | 2014-04-14 | 2015-10-15 | Kyoto University | Method to compute the barycenter of a set of histograms |
| US20150326984A1 (en) * | 2014-05-09 | 2015-11-12 | Greg Haubrich | Wireless streaming to hearing assistance devices |
| CN105279376A (zh) * | 2015-10-21 | 2016-01-27 | 武汉大学 | 一种基于gpgpu的非线性非稳态复杂信号自适应分解方法 |
| CN106446432A (zh) * | 2016-09-30 | 2017-02-22 | 云翼超算（北京）软件科技有限公司 | 一种求解材料大变形的最优输运无网格方法 |
-
2017
- 2017-12-11 CN CN201711311826.XA patent/CN108053065B/zh active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20060285114A1 (en) * | 2005-06-17 | 2006-12-21 | Canon Kabushiki Kaisha | Gas detection and photonic crystal devices design using predicted spectral responses |
| US20140192046A1 (en) * | 2013-01-07 | 2014-07-10 | Ecole Centrale Paris | Method and device for elastic registration between a two-dimensional digital image and a slice of a three-dimensional volume with overlapping content |
| US20150293884A1 (en) * | 2014-04-14 | 2015-10-15 | Kyoto University | Method to compute the barycenter of a set of histograms |
| US20150326984A1 (en) * | 2014-05-09 | 2015-11-12 | Greg Haubrich | Wireless streaming to hearing assistance devices |
| CN105279376A (zh) * | 2015-10-21 | 2016-01-27 | 武汉大学 | 一种基于gpgpu的非线性非稳态复杂信号自适应分解方法 |
| CN106446432A (zh) * | 2016-09-30 | 2017-02-22 | 云翼超算（北京）软件科技有限公司 | 一种求解材料大变形的最优输运无网格方法 |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111325337A (zh) * | 2019-08-19 | 2020-06-23 | 大连理工大学 | 一种基于自编码神经网络与最优传输的生成模型 |
| CN111008668A (zh) * | 2019-12-09 | 2020-04-14 | 上海熵熵微电子科技有限公司 | 一种数据奇异点的识别方法、装置、设备和介质 |
| CN111008668B (zh) * | 2019-12-09 | 2023-09-26 | 上海熵熵微电子科技有限公司 | 一种数据奇异点的识别方法、装置、设备和介质 |
| CN111950045A (zh) * | 2020-08-03 | 2020-11-17 | 大连理工大学 | 一种无光能损失的透镜设计方法 |
| CN113064272A (zh) * | 2021-03-04 | 2021-07-02 | 武汉大学 | 半离散最优传输下的光学自由曲面构造方法及系统 |
| CN113064272B (zh) * | 2021-03-04 | 2022-05-17 | 武汉大学 | 半离散最优传输下的光学自由曲面构造方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN108053065B (zh) | 2021-08-03 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108053065B (zh) | 一种基于gpu绘制的半离散最优传输方法及系统 | |
| JP5111638B2 (ja) | パラメトリック曲線をより小さなサブパッチに分割する装置およびその方法 | |
| US10140736B2 (en) | Graphics processing systems | |
| CN109685901B (zh) | 一种军标符号标绘中动态箭标绘制方法 | |
| WO2021203711A1 (zh) | 一种基于几何重建模型的等几何分析方法 | |
| Kowalski et al. | Smoothness driven frame field generation for hexahedral meshing | |
| US11380061B2 (en) | Method and apparatus for three-dimensional (3D) object and surface reconstruction | |
| US10540789B2 (en) | Line stylization through graphics processor unit (GPU) textures | |
| Ma et al. | Nonrigid image deformation using moving regularized least squares | |
| Yu et al. | Geodesics on point clouds | |
| CN108022005A (zh) | 一种高效的全离散最优传输方法 | |
| KR101032397B1 (ko) | 구면 좌표계를 사용하는 3차원 형상 표현장치 및 방법 | |
| CN111402422B (zh) | 三维表面重建方法、装置和电子设备 | |
| US11392806B2 (en) | Differentiable rasterizer for vector font generation and editing | |
| Jahanshahloo et al. | Reconstruction of 3D shapes with B-spline surface using diagonal approximation BFGS methods | |
| CN115713581A (zh) | 一种动态模型生成方法、装置及设备 | |
| CN113313834B (zh) | 一种基于几何可解释性的点云生成方法 | |
| CN115035233A (zh) | 全球多分辨率空间网格球面光线投射可视化方法和装置 | |
| CN115375847A (zh) | 材质恢复方法、三维模型的生成方法和模型的训练方法 | |
| CN112508007B (zh) | 基于图像分割Mask和神经渲染的空间目标6D姿态估计方法 | |
| CN114266879A (zh) | 三维数据增强、模型训练检测方法、设备及自动驾驶车辆 | |
| US11574444B1 (en) | Multi-characteristic remeshing for graphical objects | |
| Wang et al. | Robust and accurate optimal transportation map by self-adaptive sampling | |
| CN111445585B (zh) | 基于对偶凸包域的三维形状对应方法及装置 | |
| Huo et al. | A network flow approach in finding maximum likelihood estimate of high concentration regions |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |