辽宁师范大学数学学院，辽宁 大连

收稿日期：2022年11月26日；录用日期：2022年12月21日；发布日期：2022年12月30日

摘要

本文研究基于虚拟元方法求解带混合边界条件的线弹性问题。首先给定弹性力学的基本方程和二维混合弱对称形式的线弹性问题，并通过变分原理得到原问题的变分形式。其次通过定义局部刚体运动空间，构造虚拟元空间。然后对于给定的每一个变量的近似得到方法的收敛性，并通过已知的引理和不等式验证收敛性，同时给出原问题的误差估计。最后根据悬臂梁问题验证理论分析的有效性与可行性。

关键词

混合形式线弹性问题，虚拟元方法，弱对称形式，混合边界条件

Virtual Element Method for Linear Elasticity Problem

Xinbo Cheng, Xiaokun Liang, Junchi Ma^*

School of Mathematics, Liaoning Normal University, Dalian Liaoning

Received: Nov. 26^th, 2022; accepted: Dec. 21^st, 2022; published: Dec. 30^th, 2022

ABSTRACT

In this paper, the virtual element method is used to solve the linear elastic problem with mixed boundary conditions. Firstly, the basic equations of elasticity and the two dimensional mixed weakly symmetric formulation linear elasticity problems are given. The variational form of the original problem is obtained by the variational principle. Secondly, the virtual element space is obtained by defining the local rigid motion space. Then the convergence of the method is obtained for the approximation of each given variable, and the convergence is verified by the known lemmas and inequalities. At the same time, the error estimate of the original problem is given. Finally, the validity and feasibility of the theoretical analysis are verified by the cantilever beam problem.

Keywords:Mixed Formulation Linear Elasticity Problem, Virtual Element Method, Weakly Symmetric Formulation, Mixed Boundary Condition

Copyright © 2022 by author(s) and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 引言

虚拟元方法 [1] 作为有限元方法的推广，具有网格剖分的灵活性，可以使用多边形或多面体剖分，很大程度上降低了网格生成的难度，受到研究者的普遍关注。目前已被应用于求解多种问题，比如线弹性问题 [2]、膜壳问题 [3] 和Steklov特征值问题 [4] 等。本文研究应用虚拟元方法求解线弹性问题。

线弹性问题是弹性力学领域一个重要的研究课题，已有多种数值求解方法，比如，有限元方法 [5]、弱有限元方法 [6]、虚拟元方法 [7] 等。在线弹性问题中，准确地计算应力往往比计算位移更有用。在标准公式中，应力可以通过位移后处理计算，由弹性力学方程可知，这个计算需要位移的微分，从数值的角度来看，微分会使解的精确度损失，进而考虑混合形式线弹性问题，将应力张量与位移张量同时作为主要变量，但应力张量的对称性使原问题的离散化变得更加困难，为解决这一困难，可以弱化应力张量的对称性。弱化应力张量的对称性有多种方法，本文研究引入一个新的变量—旋转张量，进而对原问题进行虚拟元离散。文献 [8] 中在齐次Dirichlet边界条件下讨论混合弱对称形式线弹性问题，本文将边界条件推广为Dirichlet和Neumann混合边界条件，在该边界条件下对二维混合弱对称形式线弹性问题进行研究。

本文将使用如下符号。 $L^2\left(\mathcal{O}\right):=\left\{v:{\displaystyle {\int }_{\mathcal{O}}{\left|v\right|}^2\text{d}x}={\Vert v\Vert }_{L^2\left(\mathcal{O}\right)}^2<+\infty \right\}$ 表示定义在 $\mathcal{O}$ 上的平方可积的可测函数构成的空间， ${\left(\cdot ,\cdot \right)}_0$ 是 $L^2\left(\mathcal{O}\right)$ 上的内积， ${\Vert \cdot \Vert }_0$ 是 $L^2\left(\mathcal{O}\right)$ 上的范数。 ${\mathbb{P}}_k\left(\mathcal{O}\right)$ 为 $\mathcal{O}$ 上次数不超过k (k为非负整数)的多项式空间。 $H\left(\text{div};\mathcal{O}\right):=\left\{\tau :\tau \in {\left[L^2\left(\mathcal{O}\right)\right]}^{2\times 2}\text{,}\nabla \cdot \tau \in {\left[L^2\left(\mathcal{O}\right)\right]}^2\right\}$。定义标量积 $\sigma :\tau ={\displaystyle \underset{i,j}{\sum }{\sigma }_{ij}{\tau }_{ij}}$。定义 $as\left(\gamma \right)=\frac{1}{2}\left(\gamma -{\gamma }^{\text{T}}\right)$。

本文第二节给出弹性力学在静力学、几何学和物理学方面的基本方程，并给出二维混合形式线弹性问题，将其限制在Dirichlet和Neumann混合边界条件下。第三节给出虚拟元方法的解决，首先根据变分原理得到混合弱对称形式线弹性问题的变分形式，根据文献 [9] 证明了变分形式的解的存在唯一性。其次构造虚拟元空间，然后根据虚拟元空间得到变分形式的离散格式，最后对该方法的收敛性进行分析证明，并给出误差估计。第五节通过悬臂梁问题验证理论分析的收敛性。

2. 问题模型

弹性力学的基本方程为，在静力学方面建立平衡方程 $\nabla \cdot \sigma =f$。在几何学方面，建立几何方程

$\epsilon \left(u\right)=\frac{1}{2}\left(\nabla u+\nabla u^{\text{T}}\right)$。在物理学方面建立本构方程 $\sigma =ℂ\epsilon \left(u\right)$。其中 $\sigma$ 为应力张量， $u$ 为位移张量， $\epsilon$ 为

应变张量， $f\subset {\left[L^2\left(\Omega \right)\right]}^2$ 是外力， $ℂ$ 为正定、有界、对称的矩阵。

本文考虑的混合形式线弹性问题是线性各向同性的，找到 $\left(\sigma ,u\right)$，使得

$\{\begin{array}{ll}\nabla \cdot \sigma =f,\hfill & 在\Omega 内,\hfill \\ \sigma =ℂ\epsilon \left(u\right),\hfill & 在\Omega 内,\hfill \\ u=g_D,\hfill & 在{\Gamma }_D上,\hfill \\ \sigma \cdot n=g_N,\hfill & 在{\Gamma }_N上,\hfill \end{array}$ (1)

其中 $\Omega$ 为 ${ℝ}^2$ 上的有界凸多边形区域，设 $\Gamma$ 为 $\Omega$ 的边界，且边界 $\Gamma$ 的两个子集 ${\Gamma }_D,{\Gamma }_N$ 满足 $\Gamma ={\Gamma }_D\cup {\Gamma }_N$，${\Gamma }_D\cap {\Gamma }_N=\varnothing$。 $n$ 为边界 ${\Gamma }_N$ 的单位外法向量。 $ℂ\epsilon \left(u\right)=2\mu \epsilon \left(u\right)+\lambda tr\left(\epsilon \left(u\right)\right)I$，这里 $tr(\cdot )$ 是张量的迹，I是恒等张量，Lamé常数 $\lambda ,\mu$ 定义如下

$\lambda =\frac{E\upsilon }{\left(1+\upsilon \right)\left(1-2\upsilon \right)},\mu =\frac{E}{2\left(1+\upsilon \right)}$，

其中E为弹性模量， $\upsilon$ 为泊松比。

3. 虚拟元方法

3.1. 变分形式

构造容许函数空间 $\Sigma :=H\left(\text{div};\Omega \right)$，$U:={\left[L^2\left(\Omega \right)\right]}^2$，$X:=\left\{\gamma :\gamma \in {\left[L^2\left(\Omega \right)\right]}^{2\times 2},as\left(\gamma \right)=\gamma \right\}$。

引入一个新的变量 $\omega \in X$，根据变分原理，求得线弹性问题(1)的变分形式为：找到 $\left(\sigma ,u,\omega \right)\in \Sigma \times U\times X$，使得

$\{\begin{array}{l}a\left(\sigma ,\tau \right)+b\left(\tau ,u\right)+c\left(\tau ,\omega \right)=\langle g_D,\left(\tau \cdot n\right)\rangle ,\tau \in \left\{\tau \in \Sigma :{\tau \cdot n|}_{{\Gamma }_N}=0\right\},\hfill \\ b\left(\sigma ,v\right)=\left(f,v\right),v\in U,\hfill \\ c\left(\sigma ,\theta \right)=0,\theta \in X,\hfill \\ {\left(\sigma \cdot n\right)|}_{{\Gamma }_N}=g_N,\hfill \end{array}$ (2)

其中 $a\left(\sigma ,\tau \right):={\displaystyle {\int }_{\Omega }{ℂ}^{-1}\sigma :\tau \text{d}\Omega }$，$b\left(\sigma ,u\right):={\displaystyle {\int }_{\Omega }\left(\nabla \cdot \sigma \right)\cdot u\text{d}\Omega }$，$c\left(\sigma ,\omega \right):={\displaystyle {\int }_{\Omega }\sigma :\omega \text{d}\Omega }$。

据文献 [9] 知变分形式是适定的，则问题(2)的解存在且唯一。

3.2. 构造虚拟元空间

假设 ${\left\{{\mathcal{T}}_h\right\}}_h$ 是区域 $\Omega$ 的一个满足如下正则性条件的剖分，对任意的剖分单元 $T\in {\mathcal{T}}_h$，$h_T$ 表示每个多边形T的直径，那么令 $h:=\underset{T\in {\mathcal{T}}_h}{\sup }{h}_T$。用 ${\mathcal{E}}_h$ 表示 ${\left\{{\mathcal{T}}_h\right\}}_h$ 的边的集合。其中 $\mathcal{l}$ 为正常数，

A1. 对于所有的边 $e\in \partial T$，有 $h_e\ge \mathcal{l}{h}_T$，

A2. T相对于半径大于等于 $\mathcal{l}{h}_T$ 的球呈星形。

根据文献 [10] 构造虚拟元空间。为了构造虚拟元空间，首先分别在剖分单元 $T\in {\mathcal{T}}_h$ 上和边界 $e\in \partial T$ 上定义刚体运动空间 $RM\left(T\right)$ 和 $RM\left(e\right)$，

$RM\left(T\right):=\left\{r\left(x\right)=\alpha +\beta {\left(x-x_T\right)}^{\perp }:\alpha \in {ℝ}^2,\beta \in ℝ\right\}$，

其中若 $x=\left(x_1,x_2\right)\in {ℝ}^2$，则 $x^{\perp }={\left(x_2,-x_1\right)}^T$。

$RM\left(e\right):=\left\{\psi \left(s\right)=c{t}_e+p_1\left(s\right)n_e:c\in ℝ,p_1\left(s\right)\in {\mathbb{P}}_1\left(e\right)\right\}$，

其中 $n_e$ 是边e的单位外法向量， $t_e$ 是边e的单位切向量，方向与e一致。

下面定义局部应力空间为

${\Sigma }_h\left(T\right):=\left\{{\tau }_h\in H\left(\text{div};T\right):\exists \alpha \in {\left[H^1\left(T\right)\right]}^2,{\tau }_h=ℂ\epsilon \left(\alpha \right);{\left({\tau }_h\cdot n\right)|}_e\in RM\left(e\right),\forall e\in \partial T;\nabla \cdot {\tau }_h\in RM\left(T\right)\right\}$，

全局应力空间为

${\Sigma }_h=\left\{{\tau }_h\in \Sigma :{\tau }_{h{|}_T}\in {\Sigma }_h\left(T\right),T\in {\mathcal{T}}_h\right\}$。

定义局部位移空间为

$U_h\left(T\right):=\left\{v_h\in {\left[L^2\left(T\right)\right]}^2:v_h\in RM\left(T\right)\right\}$，

全局位移空间为

$U_h=\left\{v_h\in U:v_{h{|}_T}\in U_h\left(T\right),T\in {\mathcal{T}}_h\right\}$。

定义局部旋转张量空间为

$X_h\left({\mathcal{E}}_h\right):=\left\{{\theta }_h\in {\left[L^2\left({\mathcal{E}}_h\right)\right]}^2:{\theta }_{h{|}_e}\in RM\left(e\right),e\in {\mathcal{E}}_h\right\}$，

全局旋转张量空间为

$X_h=\left\{{\theta }_h\in X:{\theta }_{h{|}_{{\mathcal{E}}_h}}\in X_h\left({\mathcal{E}}_h\right)\right\}$。

3.3. 虚拟元离散

基于上述剖分以及空间的定义，首先将双线性形式进行离散，

$\begin{array}{l}a\left(\sigma ,\tau \right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }{a}^T\left(\sigma ,\tau \right)}={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }{\displaystyle {\int }_T{ℂ}^{-1}\sigma :\tau \text{d}T}},{\Vert \tau \Vert }_0^2={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }{\Vert \tau \Vert }_{0,T}^2},\sigma ,\tau \in \Sigma ,\\ b\left(\sigma ,u\right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }{b}^T\left(\sigma ,u\right)}={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }{\displaystyle {\int }_T\left(\nabla \cdot \sigma \right)\cdot u\text{d}T}},\sigma \in \Sigma ,u\in U,\\ c\left(\sigma ,\omega \right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }{c}^T\left(\sigma ,\omega \right)}={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }{\displaystyle {\int }_T\sigma :\omega \text{d}T}},\sigma \in \Sigma ,\omega \in X.\end{array}$

问题(2)的离散形式为，找到 $\left({\sigma }_h,u_h,{\omega }_h\right)\in {\Sigma }_h\times U_h\times X_h$，使得

$\{\begin{array}{l}{a}_h\left({\sigma }_h,{\tau }_h\right)+b\left({\tau }_h,u_h\right)+c\left({\tau }_h,{\omega }_h\right)=\langle g_D,\left({\tau }_h\cdot n\right)\rangle ,{\tau }_h\in \left\{{\tau }_h\in {\Sigma }_h:{{\tau }_h\cdot n|}_{{\Gamma }_N}=0\right\},\hfill \\ b\left({\sigma }_h,v_h\right)=\left(f_h,v_h\right),v_h\in U_h,\hfill \\ c\left({\sigma }_h,{\theta }_h\right)=0,{\theta }_h\in X_h,\hfill \\ {\left(\sigma \cdot n\right)|}_{{\Gamma }_N}=g_N,\hfill \end{array}$ (3)

其中 $\nabla \cdot {\Sigma }_h\subseteq U_h$。

根据文献 [9] 可知，离散形式是适定的，则离散形式(3)的解存在且唯一。

3.4. 收敛性分析和误差估计

首先讨论离散形式(3)的收敛性，重点在于处理双线性形式 $a$ 和 $a_h$ 之间的关系，为了证明收敛性，首先给出以下引理。

引理1 [9] 存在一个与 $\mu$ 有关的常数 $C\left(\mu \right)$，使得对所有的 ${\tau }_h\in {\Sigma }_h\cap \ker \left({\mathbb{B}}_h\right)$，有

$a_h\left({\tau }_h,{\tau }_h\right)\ge C(\mu ){\Vert {\tau }_h\Vert }_{\Sigma }^2$，

其中 $\ker \left({\mathbb{B}}_h\right)=\left\{{\tau }_h\in {\Sigma }_h:b\left({\tau }_h,v_h\right)+c\left({\tau }_h,{\theta }_h\right)=0,v_h\in U_h,{\theta }_h\in X_h\right\}$。

引理2 [9] $\ker \left(B_h\right)\equiv \left\{{\tau }_h\in {\Sigma }_h:b\left({\tau }_h,v_h\right)=0,v_h\in U_h\right\}\subseteq \ker \left(B\right)\equiv \left\{\tau \in \Sigma :\nabla \cdot \tau =0\right\}$。

引理3 [9] 存在一个常数 $C>0$，使得

$\underset{{\text{v}}_h\in U_h,{\theta }_h\in X_h}{\inf }\underset{{\tau }_h\in {\Sigma }_h}{\sup }\frac{b\left({\tau }_h,v_h\right)+c\left({\tau }_h,{\theta }_h\right)}{{\Vert {\tau }_h\Vert }_{{\Sigma }_h}\left({\Vert v_h\Vert }_U+{\Vert {\theta }_h\Vert }_X\right)}\ge C$。

定理1令 $\left(\sigma ,u,\omega \right)\in \Sigma \times U\times X$ 是变分形式(2)的解， $\left({\sigma }_h,u_h,{\omega }_h\right)\in {\Sigma }_h\times U_h\times X_h$ 是离散形式(3)的解，对于 $\sigma$ 的每一个分段多项式近似 ${\sigma }_{\pi }\in {\left[{\mathbb{P}}_k^{\text{disc}}\left(\Omega \right)\right]}^{n\times n}$，$u$ 的每一个分段多项式近似 $u_{\pi }\in U_h$，$\omega$ 的每一个分段多项式近似 ${\omega }_{\pi }\in X_h$，以及对于 $\sigma$ 的每一个近似 ${\sigma }_I\in {\Sigma }_h$，有

${\Vert \sigma -{\sigma }_h\Vert }_{\Sigma }+{\Vert u-u_{\pi }\Vert }_U+{\Vert \omega -{\omega }_h\Vert }_X\le C\left(\mu \right)\left({\Vert \sigma -{\sigma }_{\pi }\Vert }_0+{\Vert \sigma -{\sigma }_I\Vert }_{\Sigma }+{\Vert u-u_{\pi }\Vert }_U+{\Vert \omega -{\omega }_{\pi }\Vert }_X\right)$，

其中 $C\left(\mu \right)$ 为依赖于 $\mu$ 的常数。 ${\mathbb{P}}_k^{\text{disc}}\left(\Omega \right):=\left\{q_h\in L^2\left(\Omega \right):{q_h|}_T\in {\mathbb{P}}_k\left(T\right),T\in {\mathcal{T}}_h\right\}$ 为分段多项式空间。

证明：为了简便，取 $g_D=0$。

先构造 $K_h\left(f\right)=\left\{{\phi }_h\in {\Sigma }_h:b\left({\phi }_h,v_h\right)+c\left({\phi }_h,{\theta }_h\right)=\left(f,v_h\right),v_h\in U_h,{\theta }_h\in X_h\right\}$。

设 ${\tau }_h={\sigma }_h-{\phi }_h$，${\phi }_h\in K_h\left(f\right)$，则 ${\tau }_h\in \ker \left({\mathbb{B}}_h\right)$。根据引理1，引理2以及问题(2)，(3)可得

$\begin{array}{c}C\left(\mu \right){\Vert {\tau }_h\Vert }_{\Sigma }^2\le a_h\left({\tau }_h,{\tau }_h\right)\\ =a_h\left({\sigma }_h,{\tau }_h\right)-a_h\left({\phi }_h,{\tau }_h\right)\\ =a\left(\sigma ,{\tau }_h\right)-a_h\left({\phi }_h,{\tau }_h\right)+b\left({\tau }_h,u-u_h\right)+c\left({\tau }_h,\omega -{\omega }_h\right)\\ =a\left(\sigma ,{\tau }_h\right)-a_h\left({\phi }_h,{\tau }_h\right)+c\left({\tau }_h,\omega -{\omega }_h\right)\\ =a\left(\sigma -{\sigma }_{\pi },{\tau }_h\right)+a_h\left({\sigma }_{\pi }-{\phi }_h,{\tau }_h\right)+c\left({\tau }_h,\omega -{\omega }_{\pi }\right).\end{array}$

根据Cauchy-Schwarz不等式有，

$\begin{array}{c}{\Vert {\tau }_h\Vert }_{\Sigma }\le C\left(\mu \right)\left({\Vert \sigma -{\sigma }_{\pi }\Vert }_0+{\Vert {\sigma }_{\pi }-{\phi }_h\Vert }_0+{\Vert \omega -{\omega }_{\pi }\Vert }_X\right)\\ \le C\left(\mu \right)\left({\Vert \sigma -{\sigma }_{\pi }\Vert }_0+{\Vert \sigma -{\phi }_h\Vert }_0+{\Vert \omega -{\omega }_{\pi }\Vert }_X\right).\end{array}$

利用三角不等式可得，

$\begin{array}{c}{\Vert \sigma -{\sigma }_h\Vert }_{\Sigma }\le {\Vert \sigma -{\phi }_h\Vert }_{\Sigma }+{\Vert {\tau }_h\Vert }_{\Sigma }\\ \le C\left(\mu \right)\left({\Vert \sigma -{\sigma }_{\pi }\Vert }_0+{\Vert \sigma -{\phi }_h\Vert }_{\Sigma }+{\Vert \omega -{\omega }_{\pi }\Vert }_X\right).\end{array}$ (4)

根据引理3可知，存在 ${\eta }_h\in {\Sigma }_h$，满足

$b\left({\eta }_h,v_h\right)+c\left({\eta }_h,{\theta }_h\right)=b\left(\sigma -{\sigma }_I,v_h\right)+c\left(\sigma -{\sigma }_I,{\theta }_h\right),v_h\in U_h,{\theta }_h\in X_h$，

使得

${\Vert {\eta }_h\Vert }_{\Sigma }\le C{\Vert \sigma -{\sigma }_I\Vert }_{\Sigma }$。 (5)

令 ${\phi }_h={\eta }_h+{\sigma }_I$，则有

$b\left({\phi }_h,v_h\right)+c\left({\phi }_h,{\theta }_h\right)=b\left(\sigma -{\sigma }_I,v_h\right)+c\left(\sigma -{\sigma }_I,{\theta }_h\right)+b\left({\sigma }_I,v_h\right)+c\left({\sigma }_I,{\theta }_h\right)=\left(f,v_h\right)$，

其中 $v_h\in U_h,{\theta }_h\in X_h$，故可得 ${\phi }_h\in K_h\left(f\right)$。

由式(5)和三角不等式可知

${\Vert \sigma -{\phi }_h\Vert }_{\Sigma }\le C\left(\mu \right)\left({\Vert \sigma -{\sigma }_I\Vert }_{\Sigma }+{\Vert {\sigma }_I-{\phi }_h\Vert }_{\Sigma }\right)\le {\Vert \sigma -{\sigma }_I\Vert }_{\Sigma }$。 (6)

将式(6)代回式(4)可得

${\Vert \sigma -{\sigma }_h\Vert }_{\Sigma }\le C\left(\mu \right)\left({\Vert \sigma -{\sigma }_{\pi }\Vert }_0+{\Vert \sigma -{\sigma }_I\Vert }_{\Sigma }+{\Vert \omega -{\omega }_{\pi }\Vert }_X\right).$ (7)

下面分析 ${\Vert u-u_h\Vert }_U,{\Vert \omega -{\omega }_h\Vert }_X$。根据问题(2)和问题(3)可得

$\begin{array}{l}b\left({\tau }_h,u_h-u_{\pi }\right)+c\left({\tau }_h,{\omega }_h-{\omega }_{\pi }\right)\\ =a\left(\sigma -{\sigma }_{\pi },{\tau }_h\right)+a_h\left({\sigma }_{\pi }-{\sigma }_I,{\tau }_h\right)+a_h\left({\sigma }_I-{\sigma }_h,{\tau }_h\right)+b\left({\tau }_h,u-u_{\pi }\right)+c\left({\tau }_h,\omega -{\omega }_{\pi }\right)\\ \le {\Vert {\tau }_h\Vert }_{\Sigma }\left({\Vert \sigma -{\sigma }_{\pi }\Vert }_0+{\Vert {\sigma }_{\pi }-{\sigma }_I\Vert }_0+{\Vert {\sigma }_I-{\sigma }_h\Vert }_0+{\Vert u-u_{\pi }\Vert }_U+{\Vert \omega -{\omega }_{\pi }\Vert }_X\right).\end{array}$

由引理3可得

${\Vert u-u_h\Vert }_U+{\Vert \omega -{\omega }_h\Vert }_X\le C\left(\mu \right)\left({\Vert \sigma -{\sigma }_{\pi }\Vert }_0+{\Vert u-u_{\pi }\Vert }_U+{\Vert \omega -{\omega }_{\pi }\Vert }_X\right)$。 (8)

最后，联立式(7)和式(8)，结论得证。 $\square$

接下来，给出误差估计，位移、应力和旋转张量的精确解与数值解之间的误差估计如下。

定理2 [8] 问题(3)有唯一解 $\left({\sigma }_h,u_h,{\omega }_h\right)\in {\Sigma }_h\times U_h\times X_h$，则有 $L^2$ 误差估计，

${\Vert \sigma -{\sigma }_h\Vert }_0+{\Vert u-u_h\Vert }_0+{\Vert \omega -{\omega }_h\Vert }_0\le C\left(\mu \right)\left(h^{k+1}{\Vert \sigma \Vert }_{k+1}+h^{k+1}{\Vert u\Vert }_{k+1}+h^{k+1}{\Vert \omega \Vert }_{k+1}\right)$。

4. 数值实验

考虑端部受抛物线载荷作用的悬臂梁问题。在本问题中，取梁的长度为 $L=8$，高度为 $D=4$，抛物

线载荷为 $P=-1000N$，I为梁的横截面惯性矩，对于单位厚度的矩形截面 $I=\frac{D^3}{12}$。考虑本文的线弹性问

题，给出相应的材料参数，弹性模量为 $E={10}^7$，泊松比为 $\upsilon =0.3$。其精确解见文献 [11]，其中应力的精确解为

$\begin{array}{l}{\sigma }_{11}=-\frac{P\left(L-x\right)y}{I},\\ {\sigma }_{22}=0,\\ {\sigma }_{12}=\frac{P}{2I}\left(\frac{D^2}{4}-y^2\right).\end{array}$

位移的精确解为

$\begin{array}{l}{u}_1=-\frac{Py}{6EI}\left(\left(6L-3x\right)x+\left(2+\upsilon \right)y^2-\frac{3D^2}{2}\left(1+\upsilon \right)\right),\\ u_2=\frac{P}{6EI}\left(3\upsilon y^2\left(L-x\right)+\left(3L-x\right)x^2\right).\end{array}$

接下来对区域 $\Omega =\left(0,8\right)\times \left(-2,2\right)$ 进行剖分，网格剖分依照文献 [12] 完成，图1分别是网格剖分数N为100、200、400、1600的多边形剖分图。

图1. 不同网格数对应的剖分图

图2是网格剖分数N为100、200、400、1600下位移与应力的数值解图，其中 $u_{1,h}$，$u_{2,h}$ 为位移的数值解， ${\sigma }_{11,h}$ 为应力的数值解。

为了验证方法的收敛性，需要考虑不同剖分下位移和应力的相对误差。

表1给出了网格剖分数N为100、200、400、1600下的应力以及位移的相对误差，其中 $e\left(u\right)$ 表示 $L^2$ 范数下位移的相对误差，即 $e\left(u\right)={\Vert u-u_h\Vert }_0$，其中 $u$ 为位移的精确解， $u_h$ 为位移的数值解。 $e\left(\sigma \right)$ 表示 $L^2$ 范数下应力的相对误差，即 $e\left(\sigma \right)={\Vert \sigma -{\sigma }_h\Vert }_0$，其中 $\sigma$ 表示应力的精确解， ${\sigma }_h$ 表示应力的数值解。

图2. 不同网格数对应的位移与应力的数值解

表1. 不同网格数对应的位移与应力的相对误差

根据表1所给出的不同网格剖分数下的误差数据可以明显看出，随着剖分数的增加，应力以及位移的相对误差都呈现逐渐减小的趋势，则可以看出本文所讨论的方法收敛性较好，同时也验证了误差估计的有效性。

5. 总结

虚拟元方法可以看作是有限元方法对一般多边形和多面体单元的扩展。本文研究应用虚拟元方法求解具有Dirichlet和Neumann混合边界条件的线弹性问题，文中所考虑的线弹性问题是二维混合弱对称形式的线弹性问题。从给定弹性力学方程和二维连续问题出发，通过变分原理得到连续问题的变分形式，然后根据定义局部刚体运动空间构造位移、应力和旋转张量的虚拟元空间，通过构造的虚拟元空间得到变分形式的离散格式，通过引理以及三角不等式和Cauchy-Schwarz不等式证明方法的收敛性，并给出理论分析的误差估计。最后以悬臂梁问题验证方法的收敛性与误差估计的有效性为结尾。

文章引用

程新博,梁晓坤,马俊驰. 面向线弹性问题的虚拟元方法
Virtual Element Method for Linear Elasticity Problem[J]. 应用数学进展, 2022, 11(12): 9103-9111. https://doi.org/10.12677/AAM.2022.1112960

参考文献