肥皂膜及其他变分问题建模概述

2018年 9月 4日

肥皂膜、悬链线电缆和光束有什么共同点?它们都有着使某些数量最小化的行为方式,这类问题普遍存在于生物学、经济学、弹性理论、材料科学和图像处理等科学和工程领域。我们可以使用 COMSOL Multiphysics®软件中的内置物理场接口模拟很多这样的问题,在变分问题和约束系列博客中,我们将向您展示如何使用基于方程的建模功能解决变分问题。

函数的最小值

在初等微积分中,我们为单变量或多变量函数寻找最优值,即,寻找一个单一的数字或一组有限的数字\mathbf{x}使函数f(\mathbf{x})最小(或最大)。在变分法中,我们寻找一个函数u使泛函E[u(x)]最小(或最大)。从某种意义上说,我们可以认为这是无限维优化。粗略地讲,就是泛函取一个函数并返回一个数值。例如,定积分就是一个泛函。

在工程问题中,这些泛函通常代表某种能量。例如,在弹性理论中,我们可以通过最小化总势能来找到平衡解。这个术语经常被用在其他变分问题上,例如变分图像处理。我们称泛函为“能量”——即使它在物理上不代表通常意义上的能量。

Model geometry for soap film between two rings.
两个环之间的肥皂膜。

考虑 yz 平面上两个圆环之间的肥皂膜,其中心在x轴。我们希望找出这样一个函数u_x,当我们绕 x 轴旋转它时,可以获得肥皂膜的形状。这个函数实际上是下面这个泛函的最小值,

(1)

E[u(x)] = \int_a^b u(x) \sqrt{{1+u^{\prime}(x)^2}}dx.

更一般地说法是,在变分法中,我们正在寻找一个函数u(x),使泛函的值最小,

(2)

E[u(x)] = \int_a^b F(x,u,u^{\prime},u^{\prime\prime},\ldots)dx.

大多数工程问题处理的是最多包含一阶导数的泛函。在本系列博客的开始,我们将在一个空间维度上关注这些问题。稍后,我们将推广到高维度、高阶导数和几个未知数。最后,由于最大值和最小值本质是一样的,因此我们将在后续文章中只讨论最小值。

除非另有说明,否则我们将处理的泛函是

(3)

E[u(x)] = \int_a^b F(x,u,u^{\prime})dx.

求解变分问题

假设你发现自己被蒙住了眼睛进入一个山谷中,你怎么知道已经到达底部?(顺便说一句,摘掉眼罩不算)。你可能会用手和脚感觉周围的地面,如果你测试的每个区域都感觉都比你站的地方高,那么你就站在山谷的底部(至少是在一个局部的洼地)。在微积分和变分法中,同样的思想被用来检查最小值。在微积分中,我们测试相邻点,而在变分法中,我们测试相邻函数。

当且仅当对于微小值\epsilon以及每一个容许变分\hat{u(x)}都遵循以下表达式时

(4)

E[u+\epsilon \hat{u}] \ge E[u(x)]

函数u(x)使得泛函E[u(x)]最小。

不是每个变分\hat{u}是可接受的,因为每一个u+\epsilon \hat{u}都必须满足解的约束条件。例如,对于两端固定在线上的肥皂膜,我们用于与最小函数进行比较的每个函数也必须固定在线上。因此,我们只考虑那些满足\hat{u}(a) = \hat{u}(b) = 0.的变分。我们将在后面的博客文章中详细讨论约束。

假设有足够的平滑度来进行微分,方程4的必要条件是

(5)

\frac{d}{d\epsilon}\bigg|_{\epsilon=0}E[u+\epsilon \hat u] = 0,

在本系列中,我们不打算讨论移动或开放边界的问题。在这种情况下,我们可以把微分移到积分内部,并应用链式法则得到

(6)

\frac{d}{d\epsilon}\bigg|_{\epsilon=0}E[u+\epsilon \hat u] = \int_a^b [\frac{\partial F}{\partial u}\hat{u} + \frac{\partial F}{\partial u’}\hat{u^{\prime}}]dx=0.

请注意,我们只改变因变量u和它的导数,而不是空间坐标x

如果您对移动边界或界面的问题感兴趣,请查看博客文章:使用 COMSOL Multiphysics®模拟自由液面的两种方法在 COMSOL Multiphysics®中用动网格为自由液面建模

方程1所示,对于肥皂膜,我们有F(x,u,u^{\prime}) = u\sqrt{1+u^{\prime}2} \Rightarrow \frac{\partial F}{\partial u} = \sqrt{1+u{\prime}^2}, \frac{\partial F}{\partial u^{\prime}} = \frac{uu^{\prime}}{\sqrt{1+u’^2}}。因此,肥皂膜的变分问题是找到u(x),以使

(7)

\int_a^b [\sqrt{1+u’^2}\hat{u} + \frac{uu’}{\sqrt{1+u’^2}}\hat{u’}]dx=0, \forall \hat{u}.

欧拉-拉格朗日方程

在经典变分法中,我们应用分部积分从变分\hat{u}中移动空间微分到解以获得欧拉-拉格朗日方程

(8)

\frac{\partial F} {\partial u}-\frac{d} {dx}(\frac{\partial F}{\partial u^{\prime}})=0,

并利用常微分方程(ODE)方法求其解。

在高维空间中,欧拉-拉格朗日方程变成了偏微分方程。

在我们的例子中,我们不需要使用欧拉-拉格朗日方程,所以不再进一步讨论它,原因是使用有限元方法来求解变分公式。例如,在 COMSOL Multiphysics 中,如果我们使用系数形式偏微分方程或者一般形式偏微分方程接口指定欧拉-拉格朗日方程,软件会在内部指定并求解相应的变分方程,既然如此,就没有必要再浪费精力吧?正如我们将在后面看到的,变分形式还提供了很自然地考虑求解域和边界条件的方式。

在 COMSOL Multiphysics®中执行一个变分问题

为了说明 COMSOL Multiphysics 中的变分问题,我们使用了弱形式偏微分方程接口。那我们如何来对解u以及相应的试函数\hat{u}求微分呢?对于后者,我们可以使用试函数算子。例如,对于肥皂膜问题,变分公式中的被积函数在弱形式偏微分方程节点输入为sqrt(1+ux^2)*test(u)+u*ux/sqrt(1+ux^2)*test(ux),如下所示。

A screenshot of the Weak Form PDE settings in COMSOL Multiphysics.
指定一个变分问题。

我们考虑一个简单的约束,即肥皂膜固定在左右两边的金属丝环上,左环和右环的半径分别是 1 和 0.9,因此我们知道主变量u在两端的值,可以使用狄利克雷边界条件节点指定这样的边界条件。出于数值的原因,在这个特殊的问题中,我们提供了初始值 1 替代默认的初始值 0。

A screenshot of the Dirichlet Boundary Condition settings in COMSOL Multiphysics.
使用狄利克雷边界条件节点指定已知的边界值。

如果我们求解,可得到下图所示的形状。

A plot showing the profile of soap film hanging between two rings.
挂在两个垂直圆环之间的肥皂膜的轮廓。

指定一个更通用的符号形式

在上面的例子中,我们手动对F分别推导了uu^{\prime}的偏微分,而利用 COMSOL Multiphysics 的符号数学能力,我们可以避免不必要的工作和潜在的错误。

A screenshot of the Variables settings in COMSOL Multiphysics.
A screenshot of the Weak Form PDE settings when using symbolic differentiation.

用符号微分形式表达减少手动工作。

变分解与直接优化

我们还可以通过直接优化求解函数极小化问题。在这种方法中,我们不需要导出变分问题,但不利的一面是,它需要更多的计算工具。例如,在 COMSOL Multiphysics 中,直接优化需要优化模块。如果您对直接优化感兴趣,请查看我们的博客文章:如何求解两点间的最速降线

后续操作

今天,我们向您展示了如何使用 COMSOL Multiphysics 软件中的弱形式偏微分方程接口,您可以在弱形式概述这篇博文中了解更多关于弱形式的信息。

在接下来的文章中,我们将展示如何添加更复杂的约束,例如点约束、分布式约束和积分约束。本系列博客将逐步介绍更高的空间维度、更高阶的导数和更多物理场。敬请关注!

同时,您可以通过下面的按钮联系我们,了解更多关于 COMSOL Multiphysics 基于方程的建模特性。

查看变分问题和约束系列中的更多博客文章


评论 (0)

正在加载...
浏览 COMSOL 博客
Baidu
map