基于Whitney p-形式的有限元法：如何求解边系数e_i？

要从已知的电场函数$\vec{E}=(E_x,E_y,E_z)$得到Whitney边元展开系数$e_i$，核心是利用伽辽金投影将$\vec{E}$投影到Whitney边元的基函数空间中，具体步骤如下：

核心原理 #

Whitney边元的基函数${\vec{W}_i}$构成了有限元空间的一组基，展开系数$e_i$满足内积正交条件：

$\int_\Omega \vec{E} \cdot \vec{W}i \mathrm{d}\Omega = \sum_j e_j \int\Omega \vec{W}_j \cdot \vec{W}i \mathrm{d}\Omega$
其中$\Omega$是计算域，$(\cdot,\cdot)=\int\Omega \cdot \mathrm{d}\Omega$表示$L^2$内积。

分场景计算方式 #

1. 正交基特殊场景 #

如果你的Whitney基函数经过正交化处理（或网格特性天然满足正交），此时$\int_\Omega \vec{W}_j \cdot \vec{W}_i \mathrm{d}\Omega=0$（当$i\neq j$），系数可直接通过单个积分计算：

e_i = (∫_Ω E·W_i dΩ) / (∫_Ω W_i·W_i dΩ)

只需对每个边$i$，计算分子（$\vec{E}$与$\vec{W}_i$的点积积分）和分母（$\vec{W}_i$的自内积）即可。

2. 一般非正交场景 #

大多数情况下基函数非正交，需要组装线性方程组求解：

• 构造全局矩阵$A$，其中元素$A_{ij} = \int_\Omega \vec{W}_j \cdot \vec{W}_i \mathrm{d}\Omega$
• 构造右端向量$b$，其中元素$b_i = \int_\Omega \vec{E} \cdot \vec{W}_i \mathrm{d}\Omega$
• 求解线性方程组$A\vec{e} = \vec{b}$，得到系数向量$\vec{e}$

关键计算细节 #

• 单元分片积分：每个Whitney边元基函数$\vec{W}_i$仅在包含边$i$的单元（如三角形、四面体）上非零，因此积分只需遍历这些单元并累加结果。
• 基函数表达式：以3D四面体单元为例，边$i$对应的Whitney基函数为$\vec{W}i = \lambda_p \nabla\lambda_q - \lambda_q \nabla\lambda_p$，其中$\lambda_p,\lambda_q$是边$i$两个顶点的形状函数；2D三角形单元的基函数形式类似。你需要将$\vec{W}i$拆解为$(W{ix},W{iy},W_{iz})$分量，再与$\vec{E}=(E_x,E_y,E_z)$做点积计算积分。
• 数值积分实现：实际计算中用高斯积分近似，比如三角形单元选3点或6点高斯积分，四面体单元选4点或5点高斯积分。在每个积分点上计算$\vec{E}·\vec{W}_i$，乘以积分权重和单元体积/面积后累加。

代码实现提示 #

1. 遍历所有边，为每个边记录包含它的所有单元及该边在单元中的局部编号。
2. 对每个单元，根据局部边编号生成对应的Whitney基函数分量。
3. 在单元上执行数值积分，分别累加至$b_i$（右端项）和$A_{ij}$（矩阵元素，$j$为单元内的边）。
4. 利用对称正定矩阵的特性（$A$是对称正定的），选用共轭梯度法或Cholesky分解等高效求解器求解线性方程组。

内容的提问来源于stack exchange，提问作者Kristoffer Lindvall