2020年春季学期，疫情肆虐，时局动荡，纽约成为风暴中心。全球教育都转成了线上模式，我在隔离的情况下讲授了“最优传输理论和算法”的课程。在这一特殊的历史时期，每个人都得以摒除世事纷扰，沉心静气，雕琢学问，回顾人生，提炼思想。刚刚结课，回想25年来学习这一理论的心路历程，迂回崎岖，百感交集。对于丘先生的高瞻远瞩，胸襟万里，无比钦佩；对于丘先生的循循善诱，诲人不倦，由衷感激。特此撰写此文，或许对于年轻学子有所裨益。

1995年，我刚刚来到哈佛大学，当然学问处世上懵懂无知，却又磨拳檫掌，激情万丈。当时在哈佛大学学习了计算机系的课程之余，我又去旁听丘先生的拓扑课程，同时在麻省理工人工智能实验室学习Berthold Horn教授的“机器视觉”课程。Horn教授以“shape from shading”方法驰名天下，他在课堂上亲口讲授了如何用偏微分方程来描述成像过程，通过求解双曲方程来从图像反解几何。他又提出用高斯曲率来表示凸曲面的扩展高斯映射方法。这一方法其实等价于微分几何中的闵可夫斯基（Minkowski）问题。通过学习Horn教授的课程，我感到当时计算机视觉的研究方法与经典物理方法相类似，用微分方程来描述自然现象，通过解析方法求得答案。为了求解Minkowski问题，我向丘先生请教。丘先生非常兴奋地给了我他的经典论文，他在文章中证明了任意维Minkowski问题解的存在性和唯一性，并指点我去读Alexndrov和Pogorelov的经典文章和著作。多次讨论中，丘先生教给我解蒙日-安培（Monge-Ampere）方程的算法，这应该是丘先生教给我的第一个算法。

蒙日-安培方程的理论非常艰深，当时我向师兄、师姐们求助。数学系的师兄师姐都认为这是非常有深度、非常有价值的方向，值得不遗余力地努力追求。计算机系的师兄师姐们则认为，即便能够求解Minkowski问题，这种方法只能表达凸曲面，适用范围过于狭窄，因此不值得大费周章。很快，我就认识到基础科学家与工程学者的价值观念具有巨大的区别。丘先生等数学家终身追求的是大自然至善至纯的真理，美轮美奂的几何结构，而很少考虑近期人类社会的应用；应用工程学者追求的是切实可用的技术，而非长远玄虚的理论。数十年来，通过亲身实践，我既体会到自然结构带来的强烈精神享受，和超越世俗的永恒与和谐；也体会到技术进步带来的巨大震撼，和对物质强有力的掌控感。培育的博士生中，也是各有自身的价值标准，和审美倾向。很多数学出身的学生，在博士后期对实用技术出现了强烈的追求；很多从工业界回到大学的学生，对于缺乏智力深度的庸常生活难以忍受，对深刻抽象的理论产生了由衷的神往。

1995年的计算机视觉研究方式大量依赖于几何和偏微分方程方法，但是那时芒福德教授和朱松纯师兄对于Horn的研究方法提出了批评。Horn 的方法（以及当时其它的 shape-from-X 方法）中能量函数的都是人为地手工设定的，或者根据物理原理来定义的， 不能反映现实世界中物体的形状和表观的复杂性， 因此质疑这种方法的实用性。他们提出应该用统计方法来研究视觉，通过从真实图像中的统计特征， 用最大熵方法自动”学习“ 出能量函数的结构和形式， 也就是首次提出来用数据来决定模型的形式， 并且真正系统地用统计为计算机视觉的发展奠定了理论框架。当时，视觉领域几乎很少有人具有如此前瞻的想法，朱松纯师兄苦口婆心地劝我放弃微分几何的想法，去和他用统计方法做研究。当他讲到他们正在创建的统计视觉理论框架时，双眼发亮，激情难抑，“这就是一座金矿”！朱师兄的战略思想超越时代，发展过程中也是历经曲折。在机器学习处于萌芽状态的时候，在多次私下交谈里沈向洋老兄胸有成竹地预言机器学习将会大行其道，也是极具前瞻性。依随网络的发展，海量数据的积累成为现实，硬件算力迅猛提高，统计学习方法终于乘风而起。25年后，计算机视觉领域，偏微分方程的解析方法日益被边缘化，统计学习、深度学习方法几乎一统天下。

但是，历史发展的模式总是否定之否定，各个学科都在螺旋上升。近些年来，深度学习方法狂飙突进，如日中天，横扫计算机科学几乎所有领域。以前，人们为了描述一个自然物理现象，前仆后继，历经数代人努力建立了偏微分方程，无数数学家为了证明解的存在性、唯一性、正则性、稳定性，呕心沥血，苦心求索。但是在目前的深度学习框架下，我们可以将偏微分算子看成是特定函数空间中的可逆映射，即初值条件和解之间存在因果关系，那么就可以用深度神经网络去逼近这个可逆映射。这种所谓端到端的学习，省略了人类的智力加工和深刻理解这一过程，技术的进步可能带来普通人类成员智力的退化。

从基础理论层面来看，各种知识、数据规律在深度学习框架下都被表示成为概率分布，深度学习的核心目标是学习各种各样的概率分布。因此，深度学习本质上是在所有可能的概率分布组成的无穷维空间中做优化。例如传统的极大似然估计，最大后验概率，最大熵原则等等，都归结为概率空间的优化问题。那么，我们需要为这个由概率分布组成的无穷维空间建立拓扑，建立微积分。更为严谨的，假设

是底空间，它可以是欧氏空间的子集，也可以是一个黎曼流形，或者抽象的距离空间等等，其上所有可能的概率测度构成的空间为

的拓扑，如此才能谈论一族概率测度的收敛性，例如一个计算方法得到的离散概率测度如何收敛到真解；我们需要定义

中的距离，如此我们才能衡量两个概率分布之间的异同，例如GAN中判别器计算生成概率和真实概率之间的距离；我们需要定义

中的测地线，如此才能预测一个统计系统的行为，例如在通常情况下，熵会自然增大，那么所有时刻的概率分布构成了

空间中的一条轨迹，通常的热力学扩散对应着

中的一条测地线；我们需要定义绝对微分，如此我们才能在

中做优化，例如带有高斯噪声的最大似然法本质上是在

空间中极小化熵正则化的Wasserstein能量。

以上的讨论意味着，为了将深度学习奠定理论基础，我们需要从黎曼几何的角度来研究所有概率组成空间

的测地线、黎曼度量、平行移动、绝对微分等基本概念，发展概率空间的几何变分法。从这个角度来看，计算机视觉统计方法的发展又回归到无穷维黎曼几何的框架。而概率空间

的黎曼几何就是最优传输理论，而最优传输理论的核心偏方程之一就是蒙日-安培方程。丘先生25年前教会我的第一个算法就是求解蒙日-安培方程！历史再一次证明基础数学的前瞻性。丘先生曾经多次教诲我们：大自然的规律与人类经济政治无关，人类文明的真正进展是植根于对自然结构的深刻理解之上。

图0.弥勒佛计算实例，从共形参数化到保面积参数化；最后一行显示势能函数及其勒让德变换。进一步请看文后视频。识别左上图中的二维码，可以看到线上演示，并且下载数据和程序。

这学期，我为学生们讲解了最优传输理论，并且布置了实现经典算法的课后作业，其中就包括丘先生当年教给我的算法。这一算法具有强烈的几何直观，同时也具有根本的重要性。假设

是平面的凸区域，配备有两个测度

。我们希望找到一个光滑映射

，同时使得总传输代价最小，

这里代价函数为欧氏距离的平方，

这等价于求解经典的蒙日-安培方程得到凸函数（Brenier势能）：

的梯度给出最优传输映射

。最优传输映射的代价即为

的测地线由McCann平移给出：

给定测度处的切空间是所有梯度场构成的空间，两个切向量之间的黎曼内积为

由此，我们看到求解蒙日-安培方程具有根本的重要性。

我们在课上用多种理论证明蒙日-安培方程解的存在性和唯一性：

凸代价函数下的Brenier理论

不动点理论和区域不变性原理来证明Alexandrov定理解的存在性，Brunn-Minkowski 不等式证明唯一性

Perron方法证明Monge-Ampere方程解的存在性，凸函数比较定理证明解的唯一性

丘-顾-罗-孙定理【1】，用几何变分法证明解的存在唯一性

这些理论各有千秋，各有侧重和适用范围，揭示了定理的不同侧面。从计算角度而言，丘-顾-罗-孙的几何变分方法最为便捷。我们的作业就是编程实现这种方法。识别图0左上角中的二维码，可以观看线上演示，下载演示执行程序，骨架源程序和测试数据，从而可以开始动手编程。

图1. 输入是单连通曲面，带有共形参数化。

图1显示了算法的输入，一张三维单连通曲面带有一条边界

，曲面的总面积经过缩放变换，归一化成

。由曲面Ricci流方法，曲面共形映到平面圆盘

，即所谓的黎曼映照。我们在平面圆盘上放置很多小圆，拉回到曲面之后，小圆依然是小圆，但是大小发生变换。这意味着这一映射是保形的。黎曼映照将曲面的面积元前推到单位圆盘上，得到一个测度

。这一测度可以用离散的狄拉克测度之和来表示。曲面用三角网格表示，顶点映射到平面圆盘上的坐标为

。与每个顶点相邻的三角形3D面积之和的1/3记为

，这里狄拉克测度的总和等于单位圆盘的面积

图2. 输出是从曲面到单位圆盘的保面积映射。

图2显示了算法的输出。单位圆盘上的勒贝格测度（均匀分布）作为源测度

，我们的目的是求最优传输映射

下面我们用一个极度简化的模型来解释计算流程，在实际的算例中，三角网格大概有5万左右面片。图3显示了最优传输映射对应的势能函数，即蒙日-安培方程的离散解(Alexandrov解)。对应每个顶点，我们构造一个三维空间中的平面

是未知变量，满足归一化条件

。这些平面的上包络（upper envelope）构成了势能函数，

上包络向平面圆盘投影，得到圆盘的一个胞腔分解，Power Diagram，

这给出了离散的最优传输映射

，每个胞腔映射到相应的顶点。

图4. 势能函数的勒让德对偶

图4显示了势能函数的勒让德对偶，

在实际计算中，我们是先计算勒让德对偶

，因为勒让德对偶的计算归结为经典的计算几何问题：求取离散点集的凸包（convex hull）。每个平面

，势能函数的勒让德对偶等于这些对偶点的凸包：

. 目前存在大量的开源软件可以计算凸包。

如图5所示，勒让德对偶

至此，问题归结为如何计算未知高度

。我们定义高度向量的可容许空间，

根据我们的理论【1】，未知高度

向量是下面凸能量在可容许空间中的唯一全局最小点，

代表Delaunay三角剖分中连接顶点

表示在Voronoi图中的对偶边，

代表边长。同时，对角线上元素为，

由梯度向量和Hessian矩阵，我们可以用牛顿方法来优化Alexandrov能量。起始条件设置如下：初始的高度向量为：

最后一步，我们计算Power Voronoi图中每个胞腔的重心，将加权Delaunay三角剖分的每个顶点映到相应胞腔的重心，并将Power Diagram胞腔重心与加权Delaunay同样三角剖分，如此得到分片线性映射，

识别图0中的二维码，下载程序和数据包。在这个项目中，三角网格用半边结构表示（Halfedge），骨架源码包含一个半边结构的库，同时包含基于OpenGL的显示代码，以方便学生们视觉查错。凸包或者Power Voronoi图可以用CGAL或者QHull开源代码完成，骨架源码中给出了CGAL的接口。希望学生们在完成作业的时候，也显示中间计算过程，例如势能函数和勒让德对偶函数，从而加深理解，而非用最为便捷的Power Voronoi图计算。

这个项目貌似简单，实际上有很高难度：难度之一在于Power Voronoi胞腔的剪裁，Hessian矩阵元素的剪裁；难度之二在于优化过程中，高度向量

要时刻在可容许空间之中，不容走出半步。图7显示了弥勒佛的计算过程，我们看到弥勒佛头部的变化非常显著。在计算过程中，我们看到三角剖分的组合结构在不停地变化，这显示了蒙日-安培方程的强烈非线性。如图8所示，相对于人脸曲面而言，对于像手这样几何复杂的模型，计算的稳定性更加具有挑战性。这里显示的计算实例，每个顶点的相对误差都在1的负10次幂之下。

图8. 第一行手的模型，中间行共形参数化，保面积参数化，下面一行，势能函数与勒让德对偶函数。

图9. 第一行大脑皮层的模型，中间行共形参数化，保面积参数化，下面一行，势能函数与勒让德对偶函数。

目前最优传输在工程领域被广泛应用，但是由于计算的复杂性和蒙日-安培方程的强烈非线性，大量的工作中最优传输映射都是由各种方法来近似。很多优化方法也无法保证收敛到全局最优。很多时候，我们需要能够计算出ground truth，来定量比较近似方法和真解之间的差异。这次作业的方法提供了一个严密测量的工具。

理论上，这种方法对于任意维数都成立，但实际上高维三角剖分的数据结构（单纯复形）的复杂度随维数而指数增长。这为实际计算带来巨大的困难。采用粗糙的近似方法，会使得概率测度之间的距离不再精确。在计算机技术发展的历史中，这种情况屡见不鲜，最后往往是硬件技术的提升逐步解决了计算的困难。

同时，从最优传输的计算角度而言，依然存在很多开放的理论问题。例如，如何将传统的统计学习方法用最优传输理论来诠释，如何用最优传输理论来提高深度学习的可解释性，运算性能。如何分析各种最优传输映射的算法，收敛速率，误差界，稳定性等等。

下次作业，我们实现球面的最优传输映射算法；第三次作业，我们用最优传输映射来构造一个生成模型，这一模型避免模式崩溃。今年暑假，我有可能在线上讲解最优传输的基础数学理论，依然会留这些作业，欢迎有兴趣的同学参加。

深度学习本质上是在学习某种概率分布，因此深度学习主要是在所有可能的概率测度构成的空间

中进行优化。最优传输理论为

赋予了黎曼度量结构，定义了微积分，从而为概率空间中的优化奠定了理论基础。

中的绝对微分，测地线，距离等基本概念，最终归结为求解蒙日-安培方程。

25年前，丘先生教给我求解蒙日-安培方程的几何方法，这次作业就是动手实现这一基于几何变分法的算法。25年后这些理论的重要性终于被工程技术领域所认识，这再一次体现了基础科学的前瞻性，和丘先生独到深邃、超越时代的眼光。能够在求学伊始，就得到导师点化，学习到如此博大精深的理论精髓，我感到无比的荣幸，对于丘成桐先生更是充满无限感激和敬意。

在博士期间，恰逢三维技术飞速发展的时代，激光扫描、核磁共振、CT断层技术开始涌现，GPU也开始大规模普及。当时计算机科学所用的拓扑几何比较初等，基本是欧几里得几何和线性代数，无法大规模处理复杂的几何数据。时代的发展需要将微分几何、黎曼几何、共形几何、代数拓扑等现代数学引入到计算机科学当中。丘先生教了我很多数学理论，指导并亲身参与了很多研究课题。我们将曲面调和映射，度量曲面的共形结构，黎曼面单值化，全纯线丛等非常深刻的几何理论转化成实用的计算机算法，系统地创立了计算共形几何这一跨领域学科。当时计算科学家中几乎没有人从这个高度思考问题，很多离散理论（例如离散Ricci流理论、离散单值化定理）和算法（例如全纯微分，四边形网格Abel-Jacobi条件）在历史上都是我们首创。在创立计算共形几何的过程中，我们也是历尽曲折，饱受磨砺，数十年的坚韧不拔终于使得这一领域的成果广泛应用于工程和医学领域。

过去十年间，丘先生和我将我们在计算共形几何方面的研究成果加以总结，梳理汇集成《计算共形几何-理论篇》一书，由高等教育出版社和国际出版社合作出版，马上面世，欢迎感兴趣的读者关注。（识别下图中的二维码）。

图10. 《计算共形几何-理论篇》，预售二维码。

视频1. 输入曲面，弥勒佛。

视频3. 势能函数的勒让德对偶。

【2】顾险峰，丘成桐，《计算共形几何-理论篇》，高等教育出版社，International Press，2020，5.

《计算共形几何-理论篇》的二维码。

最优传输映射的线上演示，骨架源码、测试数据的二维码。包括数十个三维曲面的最优传输映射，和球面最优传输映射，以及Ricci曲率流的演示。

请长按下方二维码，选择“识别图中二维码”，即可关注。

【老顾谈几何】邀请国内国际著名纯粹数学家，应用数学家，理论物理学家和计算机科学家，讲授现代拓扑和几何的理论，算法和应用。

最近《3D游戏与计算机图形学中的数学方法》看到双向反射分布函数，开头的各种术语把爷整懵了，看完才发现竟然前面只是介绍了原理！！！

果然我还是安安分分地套公式好了。

k表示表示入射光中被漫反射出去的部分光能，1-k则表示剩余的光能，被物体吸收或成为镜面反射。可以差不多把1-k看做金属值。

D表示漫反射颜色，我们简单的用兰伯特光照模型计算即可。

重点就在于镜面反射函数 

它的镜面反射分量如下：

F 为Fresnel因子,是反射光线的强度和颜色与入射角度之间的关系函数。

D 为微平面分布函数，返回某一给定方向微平面比例。

G 为微平面的自屏蔽情况。

我们提供了视角方向V，光源方向L，中间向量H=(V+L)/2 ,通过D筛选得到了该点微平面中法向量N=H的比例。（这里的N=H仅是说明微平面的筛选条件，公式中的N指的是该点的法向量，不一定等于H）

Cook-Torrance光照模型将物体表面分解成多个微小平面，即微平面的集合,可以产生真实感很强的镜面反射。每个微平面都被当做符合电磁反射定律的理想反射平面，微平面的斜度反映了表面粗糙度。

所以我们要做的就是选择合适的F、 D、G函数，将它们组合在一起。

我找到了，它有总结FDG的各种公式，可以看一下（其中π的倒数算在D中了）。

但是这个公式感觉好复杂。。。。所以我选择使用看起来比较简单的Schlick公式：

这两个公式中的表示光线方向和法线方向一致时的折射率，你可以直接把它当做一个菲涅尔的调节值。

在公式1中，的目的就是为了得到，可能是F0比较直观吧，表示的是两个介质之间的折射率之比

当光线从法线法线照射时 L·N=1 也就是L·H=1，因为H就是L和V的中间向量，所以V·H=1.

创建一个单独的cginc文件来放置这些函数

我们使用书中的 backmann分布函数

m表示二维粗糙度向量，T表示表面的切向量，P是中间向量H到切平面的规范化投影

当mx=my时，该式等价于同向异性公式。

在顶点的切线空间计算微平面分布值会很方便，因为在切线空间中，法向量为（0,0,1），切向量为（1,0,0）公式可简化为：

光线在反射过程中会受到其余微平面的遮挡，共有两种情况（图我是偷来的(￣ε(#￣)，第三种是反射了两次的，不用管）

一是入射光线在射向目标微平面是受到路径中其余微平面的遮挡

二是反射光线在射向观察方向时受到微平面的遮挡。

取最小值作为几何衰减系数

参数的可配置范围可能不太合理，另外金属值高时，只会有高光，因为我只简单加了环境光懒得弄环境反射了