双散度椭圆方程测度边界问题:理论、正则性与数值挑战
1. 从一道“不常规”的边值问题说起在偏微分方程的理论与应用研究中我们最常接触的边值问题比如经典的泊松方程狄利克雷问题其边界条件通常是给定了边界上的函数值。然而在实际的物理建模、金融数学乃至图像处理中我们有时会遇到一种更为“奇异”的情形边界条件不是函数而是一个测度。想象一下你不是在边界上均匀地施加影响而是在边界的一个极小的、甚至零测度的点集上比如几个孤立的点集中了“能量”或“源”这种影响该如何通过方程在区域内传播并决定解的行为这正是“测度边界条件”所要刻画的核心场景。而“双散度形式椭圆方程”则进一步增加了问题的复杂性和普适性。它不再是我们熟悉的散度形式方程其最高阶项呈现出一种双重散度的结构。这种结构天然出现在许多高阶变分问题、薄板弯曲理论以及某些复合材料力学的模型中。当“双散度形式”遇上“测度边界条件”问题就从常规走向了前沿从光滑走向了奇异。这不仅仅是理论上的兴趣更是理解许多非标准物理过程数学本质的关键。本文将深入探讨这类问题的适定性解的存在唯一性、正则性解的光滑程度以及其背后深刻的泛函分析框架。2. 双散度形式椭圆方程结构、背景与动机2.1 方程的标准形式与经典散度形式对比我们首先明确什么是“双散度形式椭圆方程”。一个典型的二阶线性椭圆方程的标准散度形式为 [ -\text{div}(A(x)\nabla u) f \quad \text{in } \Omega. ] 这里未知函数 ( u ) 的梯度 ( \nabla u ) 被一个系数矩阵 ( A(x) ) 作用后再取散度。而“双散度形式”则将其提升了一个阶数常见的一种形式是 [ \text{div}(\text{div}(M(x) D^2 u)) f \quad \text{in } \Omega. ] 其中( D^2 u ) 是 ( u ) 的Hessian矩阵二阶导数矩阵( M(x) ) 是一个四阶张量系数它作用在Hessian矩阵上然后连续取两次散度。更一般地我们可以写成 [ \sum_{|\alpha||\beta|2} D^\alpha (a_{\alpha\beta}(x) D^\beta u) f. ] 这里 ( \alpha, \beta ) 是多重指标( D^\alpha ) 表示相应的偏导数。直观上看方程中未知函数 ( u ) 的最高阶导数达到了四阶但方程本身是以二阶散度的形式组织的。为什么要研究这种形式在经典弹性力学中描述薄板横向弯曲的基尔霍夫板方程在小变形假设下可简化为 ( \Delta^2 u f )双调和方程它正是双散度形式的一个特例( M ) 为单位张量。在更复杂的各向异性或非均匀材料中系数 ( M(x) ) 就不再是简单的单位张量。因此研究变系数的双散度形式方程对于理解高阶连续介质力学模型至关重要。2.2 与四阶方程及非散度形式的关系读者可能会问这本质上不就是四阶方程吗为何要强调“双散度形式”关键在于弱解的定义和所依赖的函数空间。对于完全非散度形式的四阶方程如 ( \sum_{|\alpha|4} a_\alpha(x) D^\alpha u f )要定义弱解通常要求解具有四阶弱导数这需要函数属于 ( H^4 ) 这类光滑空间对系数 ( a_\alpha ) 的光滑性要求也极高通常需要连续。而双散度形式的结构提供了巨大的便利。我们可以通过两次分部积分将方程转化为一个只涉及二阶导数的双线性形式。具体地如果试探函数 ( \phi ) 足够光滑且在边界上具有某种消失性那么形式上的计算给出 [ \int_\Omega \text{div}(\text{div}(M D^2 u)) \phi , dx \int_\Omega \sum_{i,j,k,l} M_{ijkl} \partial_{ij}^2 u , \partial_{kl}^2 \phi , dx \text{边界项}. ] 这就允许我们在仅要求 ( u ) 的二阶导数平方可积即属于 ( H^2 ) 空间的前提下定义弱解。系数 ( M ) 也仅需是有界可测的甚至可以是间断的这极大地扩展了理论的适用范围使其能处理复合材料、分层结构等实际工程问题。因此“双散度形式”不仅仅是一种书写方式它揭示了一个更弱的积分框架在这个框架下我们可以用更少的正则性假设来处理更高阶的微分算子这是现代偏微分方程理论中的一个核心思想。3. 测度边界条件当边界数据不再“温和”3.1 从函数边界到测度边界概念的跨越对于经典的狄利克雷问题我们给定 ( u|_{\partial\Omega} g )其中 ( g ) 是一个定义在边界 ( \partial\Omega ) 上的函数例如连续函数或 ( L^p ) 函数。这里的“给定”意味着在某种意义下如迹意义边界值与函数 ( g ) 相等。测度边界条件则彻底改变了这一范式。我们不再指定边界上的函数值而是指定一个作用在边界上的测度( \mu )。这个测度 ( \mu ) 可以是一个绝对连续的测度如果有密度函数则退化为函数边界也可以是一个奇异的测度比如狄拉克测度集中在单个点、一维豪斯多夫测度集中在边界的一条曲线上或者其他分形测度。问题的提法因而变为寻找函数 ( u )使得在区域 ( \Omega ) 内满足双散度形式椭圆方程并且 ( u ) 在边界上的某种“行为”由测度 ( \mu ) 所控制。但这立刻引出一个根本性问题对于一个通常属于 ( H^2(\Omega) ) 的弱解 ( u )它的迹边界值最多属于 ( H^{3/2}(\partial\Omega) ) 这类索伯列夫空间这是一个函数空间。如何用一个可能奇异的测度来匹配一个函数答案是匹配不是在逐点或几乎处处的意义上而是在一个“对偶”或“广义解”的意义上。3.2 理解测度数据的两种经典视角非常解与归一化迹在二阶椭圆方程理论中处理测度数据已有成熟的方法主要有两种视角它们都可以推广到双散度方程。第一种视角非常解Very Weak Solution。我们不再要求解 ( u ) 本身具有很高的正则性以至于其迹可以定义。相反我们降低要求通过“对偶”的方式定义解。将方程乘以一个光滑的试验函数 ( \phi )通常要求 ( \phi ) 在边界上为零即具有紧支集然后分部积分将所有导数转移到 ( \phi ) 上。这样原方程就变成了一个关于 ( u ) 的线性泛函方程。对于测度边界条件 ( \mu )我们将其视为作用在试验函数 ( \phi ) 的边界法向导数或某种高阶法向导数上的线性泛函。最终我们寻找一个函数 ( u \in L^1(\Omega) )甚至更弱的空间使得对于所有足够光滑的试验函数 ( \phi )都有某个积分等式成立该等式中包含了区域内的积分以及边界测度 ( \mu ) 作用在 ( \phi ) 的某阶法向导数上的项。这种解称为“非常解”或“分布解”。第二种视角归一化迹Normalized Trace或法向导数测度。对于二阶方程当右端项和边界数据都是测度时解可能不属于能量空间 ( H^1 )其经典的迹可能无定义。此时可以引入“法向导数测度”的概念。对于双散度方程情况更复杂但思想类似我们考虑方程在 ( H^2_0(\Omega) )边界上函数值和法向导数都为零的空间的对偶空间中的表现。边界测度 ( \mu ) 可以被解释为解在边界上的某种“广义法向弯矩”或“广义剪力”来源于板理论中的类比它作用于试验函数的相应边界量上。通过精心选择函数空间如 ( H^2 \cap H^1_0 ) 或具有更复杂边界条件的空间我们可以将测度 ( \mu ) 作为该空间对偶空间中的一个元素从而在弱形式中自然地包含边界项 ( \langle \mu, \gamma(\phi) \rangle )其中 ( \gamma ) 是某个迹算子。在实际操作中对于双散度方程由于最高涉及二阶法向导数边界测度 ( \mu ) 通常被理解为作用于试验函数 ( \phi ) 的二阶法向导数或与之相关的线性组合上。这对应于物理上在边界点或边界曲线集中施加弯矩或点力。注意选择哪种视角取决于具体的问题设置和期望的解空间。非常解理论适用范围更广对解的正则性要求最低而归一化迹方法通常与变分不等式或障碍问题结合更紧密能给出更具物理意义的解释。4. 适定性理论框架存在性、唯一性与对偶方法4.1 函数空间设置与弱形式为了严格处理双散度方程带测度边界条件的狄利克雷问题我们需要建立正确的函数空间框架。一个自然的选择是从能量空间 ( H^2(\Omega) ) 出发。标准的齐次狄利克雷边界条件对应 ( H^2_0(\Omega) )即函数及其一阶法向导数在边界上为零。但对于非齐次问题特别是测度边界条件我们需要更精细的空间。考虑空间 ( V H^2(\Omega) \cap H^1_0(\Omega) )。这个空间中的函数在边界上函数值为零但法向导数不一定为零。对于这个空间我们可以定义两个边界迹算子( \gamma_0: V \to H^{1/2}(\partial\Omega) ) 是通常的迹但在这里恒为零。( \gamma_1: V \to H^{-1/2}(\partial\Omega) ) 是法向导数迹将函数映射到其法向导数 ( \partial u/\partial \nu )。然而对于我们的测度边界条件我们可能关心的是更高阶的边界量。实际上对于双散度方程在分部积分后自然出现的边界项涉及解 ( u ) 的法向导数 ( \partial u/\partial \nu ) 和试验函数的法向导数 ( \partial \phi/\partial \nu )以及可能出现的弯矩项与 ( u ) 的二阶法向导数有关。因此一个更合适的空间可能是 ( H^2(\Omega) ) 本身并考虑其上的两个迹( \gamma_0(u) ) 和 ( \gamma_1(u) \partial u/\partial \nu )。现在假设我们的测度边界条件作用于法向导数上即我们要求 ( \partial u/\partial \nu \mu ) 在某种广义意义下成立。由于 ( \mu ) 可能是一个奇异测度它不属于 ( H^{-1/2}(\partial\Omega) )这是 ( \gamma_1 ) 的像空间的对偶而是属于其更大的对偶空间比如有界变差测度空间 ( \mathcal{M}(\partial\Omega) )或者负指数索伯列夫空间 ( W^{-2,p} )。4.2 通过对偶问题证明存在性证明这类问题弱解存在性的一个强大工具是对偶方法或Lax-Milgram定理的变体。其核心思想是将原问题转化为一个在恰当函数空间上的变分问题。定义双线性形式对于双散度方程 ( \text{div}(\text{div}(M D^2 u)) f )其自然的能量双线性形式 ( a(u, \phi) ) 定义为 [ a(u, \phi) \int_\Omega M(x) D^2 u : D^2 \phi , dx \int_\Omega \sum_{i,j,k,l} M_{ijkl}(x) \partial_{ij}^2 u , \partial_{kl}^2 \phi , dx. ] 这里要求系数张量 ( M ) 是一致椭圆的即存在 ( \lambda, \Lambda 0 ) 使得对任意对称矩阵 ( \xi )有 ( \lambda |\xi|^2 \le M_{ijkl}\xi_{ij}\xi_{kl} \le \Lambda |\xi|^2 )。处理测度边界条件将边界条件 ( \partial u/\partial \nu \mu ) 作为线性泛函纳入。我们寻找解 ( u ) 属于一个仿射空间比如 ( u \in H^2(\Omega) ) 且满足 ( \gamma_0(u) 0 )标准的狄利克雷条件但 ( \gamma_1(u) ) 是自由的。测度 ( \mu ) 定义了一个作用于法向导数迹上的线性泛函 ( L_\mu(\phi) \int_{\partial\Omega} \partial \phi/\partial \nu , d\mu )。然而直接将其作为右端项是不协调的因为 ( a(u, \phi) ) 本身已经包含了来自分部积分的边界项。构造变分形式正确的做法是考虑一个混合边值问题的弱形式。我们寻找 ( u \in H^2(\Omega) ) 满足 ( \gamma_0(u)0 )并且对于所有试验函数 ( \phi \in H^2(\Omega) ) 满足 ( \gamma_0(\phi)0 )有以下等式成立 [ a(u, \phi) \langle f, \phi \rangle_{L^2(\Omega)} \langle \mu, \gamma_1(\phi) \rangle_{\mathcal{M}, C}. ] 这里右边第一项是区域力第二项是边界测度 ( \mu ) 作用于试验函数法向导数上的对偶配对。注意试验函数空间是 ( H^2_0(\Omega) )齐次狄利克雷条件这消去了其他由分部积分产生的边界项只留下了我们关心的法向导数项。应用抽象存在性定理上述变分形式中左端 ( a(\cdot, \cdot) ) 在空间 ( { v \in H^2(\Omega): \gamma_0(v)0 } ) 上通常是强制且连续的。右端的线性泛函 ( \phi \mapsto \langle \mu, \gamma_1(\phi) \rangle ) 的连续性取决于测度 ( \mu ) 的正则性和迹算子 ( \gamma_1 ) 的连续性。如果 ( \mu ) 是一个有界拉东测度且迹算子 ( \gamma_1: H^2_0(\Omega) \to L^1(\partial\Omega) ) 是连续的实际上需要更精细的迹定理证明 ( \gamma_1 ) 到某个 ( L^p(\partial\Omega) ) 或甚至到连续函数空间的连续性那么这个线性泛函就是有界的。随后应用Lax-Milgram定理或更一般的Babuska-Lax-Milgram定理inf-sup条件即可证明弱解 ( u ) 的存在性和唯一性。实操心得在实际证明中最关键的步骤是验证由测度 ( \mu ) 诱导的线性泛函在所选试验函数空间上的有界性。这通常需要用到索伯列夫空间到边界空间的迹定理的最佳嵌入结果。对于奇异测度如狄拉克测度需要证明试验函数的法向导数在测度支撑点上的取值是良定义的这要求试验函数具有足够的正则性例如连续因此可能需要将试验函数空间限制在更光滑的子空间如 ( C^\infty ) 紧支集函数然后通过稠密性进行延拓。5. 解的正则性分析从奇异性到局部平滑5.1 测度数据导致的必然奇异性一个必须接受的现实是当边界数据是奇异测度如点质量时解在区域内不可能处处光滑。以最简单的例子说明在单位圆盘上考虑拉普拉斯方程 ( \Delta u 0 )边界条件为在边界某点 ( P ) 处放置一个单位狄拉克测度。这个问题的解在某种广义意义下等价于该点的泊松核它在 ( P ) 点具有奇异性。对于双散度方程情况类似但更复杂奇异性可能更强。解的正则性主要取决于两个因素1)测度 ( \mu ) 的奇异性2)方程的阶数。对于二阶椭圆方程如果右端项或边界数据是狄拉克测度解通常属于 ( W^{1,p} ) 对于所有 ( p n/(n-1) )在n维空间中但不属于 ( W^{1, n/(n-1)} )更不属于 ( H^1 )。这意味着解的一阶导数在测度支撑点附近是可积的但积分幂次有限。对于四阶的双散度方程由于方程阶数更高解对奇异性有更强的“平滑效应”。直观上方程像是一个“更硬的弹簧”能将奇点的能量更有效地扩散开。然而边界上的点测度仍然会导致解在边界点附近失去正则性。通常解 ( u ) 可能不再属于 ( H^2(\Omega) )能量空间而只能属于 ( W^{2,p}(\Omega) ) 对于某个 ( p 2 )或者属于更弱的 ( W^{2,1} ) 空间。其导数可能在边界奇点附近爆炸。5.2 Calderón-Zygmund型估计与局部正则性尽管整体正则性受损但在远离奇异测度支撑的区域解仍然可以表现出很好的光滑性。这是椭圆方程局部正则性理论的体现。对于双散度方程如果系数 ( M(x) ) 足够光滑例如一致连续那么对于任何内部子区域 ( \Omega \subset \subset \Omega )即远离边界解 ( u ) 是光滑的( C^\infty )。即使系数只是有界可测通过Gehring引理或Meyers估计我们也能得到解的二阶导数在局部属于某个 ( L^{p}_{loc} ) 空间其中 ( p 2 )这比全局的正则性要好。对于边界附近但远离测度奇点的情况情况类似。如果边界是光滑的且测度 ( \mu ) 在边界的一个闭子集 ( K ) 上是奇异的那么在边界上任意一个不包含 ( K ) 的开子集上解 ( u ) 的法向导数甚至更高阶的边界迹可能是 Hölder 连续的。这种“局部平滑”现象在实际计算中非常重要它意味着数值误差或物理扰动在远离奇点的地方不会剧烈放大。分析此类正则性的标准工具是差商法和Campanato空间理论。基本步骤是对弱形式进行平移得到关于差商 ( \Delta_h u ) 的方程。利用方程的椭圆性对差商建立Caccioppoli型不等式能量估计。通过迭代或嵌入定理将这种能量估计提升为高阶可积性或Hölder连续性估计。对于测度数据问题这个过程需要在避开奇点支撑的区域内进行。一个实用的技巧是使用截断函数选取一个光滑的截断函数 ( \eta )其在奇点支撑的一个邻域外恒为1在该邻域内光滑地衰减到0。然后对 ( \eta u ) 应用正则性估计。由于 ( \eta ) 的引入会在方程中产生低阶项但这些项是光滑系数乘以 ( u ) 的低阶导数可以通过吸收技巧和已知的先验估计例如解的整体 ( L^1 ) 或 ( W^{1,1} ) 估计来控制。注意事项在进行局部正则性估计时常数通常依赖于到奇点支撑的距离。距离越近常数可能越大甚至趋于无穷。这意味着在非常靠近奇点的区域解的振荡可能无法被标准的内估计所控制。在数值计算中这提示我们需要在奇点附近进行网格加密或使用自适应方法。6. 数值逼近的挑战与有限元方法策略6.1 直接离散化的困难将测度边界条件直接纳入有限元框架并非易事。主要的困难在于泛函的离散化线性泛函 ( L_\mu(\phi) \int_{\partial\Omega} \partial \phi/\partial \nu , d\mu ) 如何用有限维的试探函数来近似如果 ( \mu ) 是狄拉克测度 ( \delta_{x_0} )那么 ( L_{\delta_{x_0}}(\phi) \partial \phi/\partial \nu (x_0) )。这要求有限元空间中的函数在点 ( x_0 ) 处的法向导数有定义。对于标准的拉格朗日有限元如 ( P1 ) 或 ( Q1 )其法向导数在节点间不连续在节点处的值依赖于网格并非一个内在定义的量。解的低正则性由于解 ( u ) 不属于能量空间 ( H^2 )伽辽金投影的最佳逼近误差估计 ( \inf_{v_h \in V_h} | u - v_h |_{H^2} ) 可能没有意义甚至无穷大。标准的先验误差分析框架失效。边界层效应在奇点 ( x_0 ) 附近解的高阶导数会急剧变化形成边界层。均匀网格会在此处产生巨大误差。6.2 正则化策略与对偶混合方法针对这些挑战有两种主流的数值处理策略。策略一测度正则化。这是最直观的方法。将奇异测度 ( \mu ) 用一个光滑函数序列 ( \mu_\epsilon ) 来逼近例如用磨光核卷积狄拉克测度得到一个在 ( x_0 ) 附近小邻域内光滑分布的“力”。然后求解正则化后的问题 [ \text{div}(\text{div}(M D^2 u_\epsilon)) f, \quad \partial u_\epsilon/\partial \nu \mu_\epsilon \text{ on } \partial\Omega. ] 对于每个固定的 ( \epsilon 0 )( u_\epsilon ) 具有更好的正则性可以用标准有限元方法如 ( C^1 ) 连续的元如埃尔米特元或Argyris元进行离散。关键点在于分析误差 ( | u - u_\epsilon | ) 以及离散误差 ( | u_\epsilon - u_{\epsilon, h} | )并选择合适的 ( \epsilon ) 与网格尺寸 ( h ) 的关系。通常需要 ( \epsilon ) 趋于零的速度慢于 ( h )以保证整体误差收敛。策略二对偶混合有限元法。这种方法避免直接处理测度边界条件。我们引入一个新的变量 ( \sigma M D^2 u )可以理解为弯矩张量将四阶方程降阶为一个二阶方程组 [ \begin{cases} \sigma - M D^2 u 0 \text{in } \Omega, \ \text{div}(\text{div } \sigma) f \text{in } \Omega. \end{cases} ] 然后对其应用混合变分形式。边界条件 ( \partial u/\partial \nu \mu ) 可以转化为对偶变量 ( \sigma ) 的某种边界积分条件。在这种形式下测度 ( \mu ) 出现在线性泛函中作用于试探函数的边界项上。离散时我们可以分别对 ( u ) 和 ( \sigma ) 选择有限元空间通常 ( u ) 用 ( H^1 ) 元( \sigma ) 用 ( H(\text{div}, \text{div}) ) 类型的元。这种方法的好处是降低了对解 ( u ) 本身正则性的要求并且能更自然地处理边界积分。然而它需要满足严格的 inf-sup 条件以确保离散格式的稳定性这对有限元空间对的选择提出了挑战。6.3 自适应网格细化无论采用哪种策略在奇点附近进行自适应网格细化都是至关重要的。基于后验误差估计子如残差型估计子算法可以自动识别解梯度或曲率变化剧烈的区域即奇点附近并在这些地方加密网格。对于双散度方程一个有效的后验误差估计子可能包含以下几项单元残差( \eta_K h_K^2 | f \text{div}(\text{div}(M D^2 u_h)) |_{L^2(K)} )对于内部单元。法向跳跃残差对于 ( C^0 ) 连续的有限元如混合法中的 ( u_h )其弯矩 ( \sigma_h M D^2 u_h ) 的法向分量在单元边界上可能不连续。这个跳跃的大小 ( \eta_E h_E^{3/2} | [\sigma_h \cdot n] |_{L^2(E)} ) 是误差的重要指示器。边界残差对于测度边界条件如果采用正则化边界残差 ( \eta_{\partial\Omega} | \mu_\epsilon - \partial u_h/\partial \nu |_{L^2(\Gamma)} ) 在奇点支撑附近会很大驱动局部细化。通过迭代求解→估计误差→标记并细化网格的过程可以在奇点附近获得指数级的收敛速度从而用相对较少的自由度获得高精度的解。实操心得在实现正则化方法时正则化参数 ( \epsilon ) 的选择需要与网格尺寸 ( h ) 协调。一个经验法则是让 ( \epsilon ) 与局部网格尺寸 ( h_{local} ) 同阶或略大。如果 ( \epsilon ) 太小而网格太粗正则化的力无法被网格有效分辨会导致数值解振荡如果 ( \epsilon ) 太大则正则化误差占主导。通常的做法是在自适应细化过程中让 ( \epsilon ) 随着网格加密而逐步减小。7. 物理应用举例集中弯矩作用的弹性薄板让我们用一个具体的物理模型来贯穿上述抽象理论承受集中弯矩作用的弹性薄板。考虑一个由均匀各向同性材料制成的薄板占据平面区域 ( \Omega )。其在外力作用下的横向位移 ( w(x, y) ) 由基尔霍夫板方程控制 [ D \Delta^2 w q(x, y). ] 其中 ( D ) 是板的弯曲刚度( q ) 是横向分布载荷。这是一个双散度形式方程的特例( M ) 为常数张量。现在假设在板边界 ( \partial\Omega ) 上的某一点 ( P ) 处施加了一个集中弯矩 ( M_0 )单位是力×长度。在经典的连续力学中边界条件通常表述为弯矩或剪力的分布。一个集中弯矩在数学上无法用函数来描述它本质上是一个狄拉克测度。因此边界条件应写为 [ \frac{\partial w}{\partial \nu} M_0 \delta_P \quad \text{on } \partial\Omega. ] 这里( \partial w/\partial \nu ) 是位移的法向导数在板理论中与转角相关。( \delta_P ) 是边界点 ( P ) 处的狄拉克测度。这个模型完美契合了我们讨论的框架方程双散度形式双调和算子。边界条件测度边界条件狄拉克测度。物理解释测度 ( \mu M_0 \delta_P ) 代表了在点 ( P ) 处集中作用的广义力弯矩。理论分析根据前面的适定性理论只要区域 ( \Omega ) 是 Lipschitz 区域且 ( M_0 ) 是有限常数这个问题在 ( H^2(\Omega) \cap H^1_0(\Omega) ) 的某种对偶意义下存在唯一的弱解或非常解。解 ( w ) 在点 ( P ) 附近具有奇异性其曲率二阶导数在 ( P ) 点趋于无穷。但在远离 ( P ) 的板内部位移 ( w ) 是光滑的。数值模拟直接使用有限元软件模拟此问题会遇到困难。标准的板单元如DKQ、DKT或基于Kirchhoff理论的 ( C^1 ) 连续元无法直接施加点弯矩。实践中工程师通常采用两种近似正则化将点弯矩 ( M_0 ) 等效为作用在 ( P ) 点附近一个微小线段 ( \Gamma_\epsilon ) 上的均布弯矩 ( M_0 / |\Gamma_\epsilon| )。这对应于用特征函数测度 ( \frac{M_0}{|\Gamma_\epsilon|} \chi_{\Gamma_\epsilon} ) 来逼近 ( \delta_P )。对偶混合法引入弯矩 ( \mathbf{M} ) 和剪力 ( \mathbf{Q} ) 作为独立变量建立混合变分形式。集中弯矩 ( M_0 \delta_P ) 自然出现在边界虚功项 ( \int_{\partial\Omega} M_0 \frac{\partial \delta w}{\partial \nu} , d\delta_P M_0 \frac{\partial \delta w}{\partial \nu}(P) ) 中。在离散时需要确保试探函数 ( \delta w ) 在点 ( P ) 处的法向导数有定义这要求网格节点必须包含 ( P ) 点并使用至少 ( C^1 ) 连续的形函数如埃尔米特元。结果解读数值解会显示在 ( P ) 点附近板的曲率非常大应力集中。位移场 ( w ) 在 ( P ) 点处是连续的但其梯度转角在 ( P ) 点可能有一个“跳跃”或急剧变化。通过后处理计算出的弯矩场 ( \mathbf{M} ) 在 ( P ) 点附近会出现极高的值这符合集中载荷下应力奇异的物理预期。自适应网格细化会显著改善 ( P ) 点附近解的精度。