COMSOL Multiphysics弱形式入门

COMSOL Multiphysics弱形式入门物理问题的描述方式有三种：1、偏微分方程2、能量最小化形式3、弱形式本文希望通过比较浅显的方式来讲解弱形式，使用户更有信心通过COMSOL Multiphysics的弱形式用户界面来求解更多更复杂的问题。

COMSOL Multiphysics是唯一的直接使用弱形式来求解问题的软件，通过理解弱形式也能更进一步的理解有限元方法（FEM）以及了解COMSOL Multiphysics的实现方法。

本文假定读者没有太多的时间去研究数学细节，但是却想将弱形式快速的应用到实际工程中去。

另外，本文也会帮助理解COMSOL Multiphysics文档中常用的到一些术语和标注方法，相关理论可以参考Zienkiewicz[1]，Hughes[2]，以及Johnson [3]等。

为什么必须要理解PDE方程的弱形式？一般情况下，PDE方程都已经内置在COMSOL Multiphysics的各个模块当中，这种情况下，没有必要去了解PDE方程和及其相关的弱形式。

有时候可能问题是没有办法用COMSOL Multiphysics内置模块来求解的，这个时候可以使用经典PDE模版。

但是，有时候可能经典PDE模版也不包括要求解的问题，这个时候就只能使用弱形式了（虽然这种情况是极少数的）。

掌握弱形式可以使你的水平超过一般的COMSOL Multiphysics用户，让你更容易去理解模型库中利用弱形式做的算例。

另一个原因就是弱形式有时候描述问题比PDE方程紧凑的多。

还有，如果你是一个教授去教有限元分析方法，可以帮助学生们直接利用弱形式来更深入的了解有限元。

最后，你对有限元方法了解的越多，对于COMSOL Multiphysics中的一些求解器的高级设置就懂得更多。

一个重要的事实是：在所有的应用模式和PDE模式求解的时候，COMSOL Multiphysics 都是先将方程式系统转为了弱形式，然后进行求解。

PDE问题常常具有最小能量问题的等效形式，这让人有一种直觉，那就是PDE方程都可以有相应的弱形式。

实际上这些PDE方程和能量最小值问题只是同一个物理方程的两种不同表达形式罢了，同样，弱形式（几乎）是同一个物理方程的第三个等效形式。

这三种形式的区别虽然不大，但绝对是很关键的。

我们必须记住，这三种形式只是求解同一个问题的三种不同形式――用数学方法求解真实世界的物理现象。

根据不同的需求，这三种方式又有各自不同的优点。

PDE形式在各种书籍中比较常见，而且一般都提供了PDE方程的解法。

能量法一般见于结构分析的文献中，采用弹性势能最小化形式求解问题是相当自然的一件事。

当我们的研究范围超出了标准有限元应用领域，比如传热和结构，这个时候弱形式是不可避免的。

化工中的传质问题和流体中的N-S方程都是没有办法用最小能量原理表述出来的。

本文后面还有很多这样的例子。

PDE方程是带有偏微分算子的方程，而能量方程是以积分形式表达的。

积分形式的好处就是特别适合于有限元方法，而且不用担心积分变量的不连续，这在偏微分方程中比较普遍。

弱形式也是积分形式，拥有和积分形式同样的优点，但是他对积分变量的连续性要求更低，可以看作是能量最小化形式的更一般形式。

最重要的是，弱形式非常适合求解非线性的多物理场问题，这就是COMSOL Multiphysics的重点了。

小结：为了理解PDE方程的弱形式，我们必须跳开常规的偏微分形式，对于积分形式要好好研究。

由于最小于能原理对比弱形式来说好理解的多，所以我们将从线弹性开始学习，依次到热传导，电流传导等问题。

这几种物理问题都有相关的能量和功率可以进行最小化。

我们将只涉及到静态问题，重点是在结构分析和更特殊的线弹性分析。

弹性静力学PDE 及其弹性能量方程在静力结构分析问题中，我们需要求解的是Navier 方程σ-∇⋅=F其中σ是应力张量，F 是体力，比如重力等。

如果不习惯用张量的形式，你也可以将张量展开写成矩阵形式。

这个方程表示了力（或者等效力）的平衡，实际上是三个方程的合并形式——3D 中每个坐标方向有一个方程。

计算区域记为Ω，其边界记为Ω∂。

应力张量σ和应变张量ε之间的关系称为本构关系，线弹性本构一般遵循胡克HOOK 定律εσεc =其中εc 是弹性张量，这个关系式说明材料的行为实际上和弹簧差不多（前提是线弹性）。

最后，我们可以将应变矢量和位移的关系表述出来ε=∇u这里u 指的是位移矢量u =(u,v,w)，其定义就是变形体上的材料点和未变形时候的位移差。

总结以上所有的方程，我们得到了一个二阶PDE 方程（Navier 方程），()c -∇⋅∇=u F需要一个边界条件来求解，()n c ⋅∇=u P其中n 是表面Ω∂的法矢，P 是边界上的面力或牵引力。

后面会介绍更多边界条件。

这个PDE 方程的弱形式为，其中v=),,(z y x v v v 称为试函数。

注意，尽管Navier 方程是一个矢量表达式，但是上面的表达式是一个标量形式。

下面介绍如何去推导以及理解弱形式。

弹性势能在结构分析中，PDE 方程及其弱形式的表达式都不太常见，相反，能量最小化形式因为其直观的表达形式用的较多。

这类问题的能量积分形式对应于总势能的最小化，即对象中存储的弹性能。

总弹性能是一个标量，可以写成：弹性能表达式同样适用于非线性问题。

在这些表达式中，我们假设体力F 为零，并忽略了边界效应。

这些影响可以在以后引入。

积分的意义是每个体积微元的内能总和，其中应力张量单位是Pa ，微元体上的应变εd 没有单位，dV 单位是体积，因此积分出来的单位应该是N·m 。

如果问题是线弹性的，则可以显式的写为：利用下面的通用公式：用应变张量ε替换上式中的标量变量x ，弹性张量εc 替换上标量常量a 。

联立上面的式子得到：我们用c 代替εc 来配合COMSOL Multiphysics 手册中的标记方式。

再提醒一次，如果你不习惯用张量，可以将张量看成是一个3×3的矩阵，点乘是一种张量的运算符号，弹性张量εc 是一个4阶张量（看上去就像4维矩阵）。

更多的标记方法可以参考COMSOL Multiphysics 的Anisotropic Structural Analysis 中的Matrix Notation 。

弹性能积分形式下的单位说明：32][m N [c]无][m d Pascal =Ω===⎰Ω单位ε 最终给出总的积分单位是N·m ――能量。

E W 的表达式就是我们通常说的能量泛函，即位移矢量u （或实际上是u 的梯度）的泛函。

这种函数的函数，而不是坐标的函数，通常被称为泛函，比单元微积分和多元微积分更加抽象。

与积分类似，我们可以说E W 就是函数u 的泛函：这好比是一个2D 的变量x ，y 的二元函数：其中x ＝),(y x ，[]a b A c d=，[]d b e =。

采用这样的类比是因为在后面我们会看到矩阵A 与有限元的刚度矩阵比较类似。

我们要说明一下函数和泛函的一些区别，古典分析中的函数概念是指两个数集之间所建立的一种对应关系，现代数学的发展却是要求建立两个任意集合之间的某种对应关系。

函数概念被赋予了更为一般的意义，通俗解释泛函指的就是“函数的函数”。

在这里定义域为Ω，泛函可以在整个定义域内进行微分积分等操作。

泛函的变量是函数，这个函数也是有容许空间的。

如果函数u 可以变化，可能会产生一些不符合物理规则的一些现象，例如结构的刚性位移等。

比如一个对u 的基本约束就是材料不能穿越本身。

在有限元分析中，泛函一般是某种能量积分，比如弹性能。

对于其他的物理场，可能是其他的能量积分，或者是一种等效于能量的标量也可以。

至于积分区域，一般由分析对象的CAD 几何区域所确定。

静态电流传导和能量的生成在静态导电问题中，PDE 方程由最基本的保守形式开始：0∇⋅=J其中J 是电流密度。

材料（或本构）模型采用欧姆Ohm 定律：r σ=B E其中E 是电场，r σ是电导率。

另外，已知：V =-∇E其中V 是静电势，综合以上式子得到()0r V σ-∇⋅∇=在COMSOL Multiphysics 中，这就是所谓的Conductive Media DC 方程。

电阻产生的热能稳态电流的能量问题是在电导体中的电阻热其中J 表示电流强度，E 代表电场强度，εσ是一个二阶电导张量(3×3)。

如果导体是金属，电导张量一般是一个对角矩阵，如果是晶体，情况就复杂多了。

尽量减少电阻产生的热量，也就是减少热损耗，是我们要研究的一个最小值问题。

如果问题是线性，则积分可以显式地写成：因为V =-∇E ，其中V 是电势，可以得到：将这个式子与结构力学中的式子进行对比，发现他们非常相似。

V 的梯度对应于位移梯度，电导率张量εσ对应于弹性张量c 。

在稳态电流和结构力学的计算过程中，张量形式都可以改写为矩阵形式。

传热PDE 方程和能量形式对于稳态传热问题，PDE 形式为：其中T 是温度，k 是热传导系数，Q 是空间分布的热源。

热能基于传热方程的典型泛函为：其中T 是温度，k 是热传导系数张量（3×3）。

泛函极小值泛函极值的概念借用了微积分中的不少方法。

本节首先会介绍函数微积分的求极值方法，接下来，我们会借用有限元中常用的术语和标注方法来推导我们熟悉的结果。

这个过程可以被看作是微积分方法的一种推广。

考虑一个多元微积分函数f，我们要求最小值：寻找x使得f(x) 最小化这里x是一个矢量，或者点的坐标。

通过微积分我们知道，这个时候首先必须求函数f的梯度。

将梯度的设置为0，我们可得到一个非线性方程组。

求解方程，我们可以得到一系列的坐标点x，如果在其中某点处的二阶倒数(一般称为Hessian矩阵)为正(或者说有正的特征值)，就说这点就是我们要求的极小点，就好像该点是整个函数的一个谷底一样。

利用Taylor展开的观点，假设已知一个最小值x，我们可以在上面施加一个小的扰动，由Taylor展开可得：这里H就是前面所说的Hessian矩阵。

现在我们用其他的方法来说明函数f在x最小。

首先，假设x是一个极值点，当添加了一个xδ后，f对于其一阶值不改变。

换句话说，如果我们在x上添加一个xδ来扰动f，其一阶Taylor级数应该为0。

这个条件应该对每个方向都是成立的，否则该点就不是极值点了。

如果上式第二项为0：对于任意小的xδ都成立，也就是：我们这里只是用一个稍微有点不同的方法得到了一个同样的结果。

但是，这只是给了我们一个极值点的信息，如果要确定其是最小极值点，必须保证第三项（二阶项）对于任意xδ都为正：只有当H的特征值都为正时，上式成立（参考线性代数）。

有可能会遇到二阶项也总为0，这个时候我们必须借助更高阶项来判断极值点。