fealpy
fealpy copied to clipboard
weihuayi
LagrangeFiniteElementSpace 中 grad_basis 不支持 TD-1 维的积分点
需要求 grad_basis() 在边(面)上的值.
能否将 grad_basis() 中的 TD=self.TD 改为类似于 basis() 中的 TD = bc.shape[-1] - 1 ?
应该不可以，关键你说的 grad_basis 在边（面）上的值的具体数学含义，是单元的梯度值在边界上的限制呢？还是边界上切身梯度值？
或者你需求背景是什么？
用速度校正格式求 N-S 方程, 中间过程有一步需要得到 <fun^n, \nabla q>_\Gamma, 其中 fun^n 是上一时间层的值 (不用管), 主要 q 是压力空间中的 test-function, \Gamma 可以简单当做 Dirichlet 边界.
所以你的 \nabla q 是单元值在 \Gamma 上的限制，对吗? 你会给出边上的积分点，然后计算 \nabla q 在些积分点上的值？
是的, face_basis() 和 basis() 两个代码是一样的, 所以我最开始的直观理解是 grad_basis() 也可以根据 bc.shape 来计算边(面) 上的梯度值.
@yoczhang 你看一下 这个代码 edge_grad_basis
可以试下
另外，你的数值离散必须要在三角形网格上做吗？
意思是 edge_basis 和 edge_grad_basis 只适用于 2d 中的边?
还有点困惑的是, fealpy 的 FEM 中 Lagrange 基函数的实现, 比如 2d 中我限制在某个边上求 basis(bc), 此时 bc 为 1d 的 Gauss-Legendre-Quadrature, basis() 也退化为 1d 的情况, 但是此时所返回的 basis(bc) 就是 2d 中在此边上的限制吗?
edge_grad_basis 中 lidx 表示什么?
我目前测试中出现错误:
gphi = np.einsum('k...ij, kjm->k...im', R, Dlambda[index, :, :]), index 中的指标超出了 Dlambda 的维数
index 和 lidx 长度应该一样， 你要计算一个单元的基函数在某条边上积分点处的函数值，你需要知道单元的编号 index 和 这条边在这个单元上的局部编号。
我问 “另外，你的数值离散必须要在三角形网格上做吗？" 的意思是： 如果你的的离散方法适用于多边形网格的，最好不用 Lagrange 有限元空间。