ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

PyNite DKMQ板单元揭秘:四边形板有限元公式推导详解

PyNite DKMQ板单元揭秘:四边形板有限元公式推导详解 PyNite DKMQ板单元揭秘四边形板有限元公式推导详解【免费下载链接】PyNiteA 3D structural engineering finite element library for Python.项目地址: https://gitcode.com/gh_mirrors/py/PyNitePyNite 是一个用 Python 编写的 3D 结构工程有限元库本文带你完整看懂其 DKMQ 四边形板单元的有限元公式推导过程从自由度布置、双线性形状函数、离散 Kirchhoff 约束到弯曲刚度矩阵的组装每一步都讲透新手也能轻松跟上 一、为什么 DKMQ 板单元值得深究PyNite 提供两种板单元它们定位不同详见文档 docs/source/plate.rst 的说明单元类型几何要求核心思想适用场景Rect矩形板必须为矩形12 项多项式弯曲函数矩形网格、快速建模Quad四边形板任意四边形DKMQ 等参数化公式厚薄板通吃、扭曲网格DKMQDiscrete Kirchhoff Mindlin Quadrilateral公式的精髓在于把 Kirchhoff 板理论的高精度弯曲行为与 Mindlin 理论的横向剪切变形嫁接在一起。这就让它既不失薄板的精度又不会在厚板时出现剪切锁定——这也是 PyNite 官方示例中评价它对厚板和薄板都能给出很准确结果的原因。实现代码位于 Pynite/Quad3D.py文件头部列出了 4 篇经典参考文献Katili 的 DKMQ/DSQ/MITC4 对比研究、Bathe、Logan、Gallagher是学习四边形板有限元公式推导的宝藏起点。二、自由度与局部坐标系板单元的关节 一个 DKMQ 四边形板单元共有16 个自由度4 个角节点每个节点带 3 个弯曲自由度——横向位移w和绕两条板面轴的转角βx、βy合计 12 个4 个边中点节点每个边中点带 1 个绕边法线的转角Δβs合计 4 个。边中点转角的存在正是 DKMQ 区别于普通双线性单元的关键它让板面斜率可以表达出二次变化从而显著改善弯曲精度。在结构组装层面每个角节点有 6 个自由度3 平移 3 转角4 个节点共 24 个单元刚度矩阵为 24×24 阶。其中绕板面法线的转动即钻孔自由度在纯弯曲理论中是无约束的PyNite 通过弱旋转弹簧刚度取其他转动刚度的 1/1000来保证数值稳定这一做法在 Pynite/Quad3D.py 的类说明和 Pynite/Plate3D.py 的ke_b方法中都有体现。单元的局部坐标系由节点顺序决定i → j方向为局部 x 轴法向量由叉积确定 z 轴。理解这套局部坐标是读懂所有板单元推导的前提PyNite有限元成员局部坐标系与截面内力方向定义三、形状函数双线性角点 不完整二次边中点DKMQ 在自然坐标系(ξ, η) ∈ [-1, 1]上插值横向位移w形状函数分两类对应 Pynite/Quad3D.py 中的N_i与P_k方法角点双线性函数i 1~4$$N_i \frac{1}{4}(1 \pm \xi)(1 \pm \eta)$$边中点不完整二次函数k 5~8$$P_k \frac{1}{2}(1 - \xi^2)(1 \mp \eta),\quad P_k \frac{1}{2}(1 \pm \xi)(1 - \eta^2)$$注意不完整二字二次项里缺少ξη交叉项。这不是偷懒而是刻意设计——去掉交叉项后边中点转角Δβs与角点转角在单元内部通过离散 Kirchhoff 约束保持协调在单元内的高斯积分点上转角与位移之间的 Kirchhoff 条件转角 位移斜率被逐点强制成立而无需像 C¹ 连续单元那样要求跨节点斜率连续。这一机制的完整符号推导可以在Derivations/DMKQ Quad Element.ipynbJupyter Notebook需用 Sympy 逐格运行中逐步查看。四、弯曲刚度矩阵的推导链条 弯曲刚度的推导是 DKMQ 公式的心脏核心链条如下对应 Pynite/Quad3D.py 中的各矩阵方法1. 剪切耦合系数 φkphi_k方法$$\phi_k \frac{2}{\kappa(1-\nu)}\left(\frac{t}{L_k}\right)^2,\quad \kappa \frac{5}{6}$$其中Lk为第 k 条边的长度t为板厚ν为泊松比。板越薄φk 越小约束越接近刚性 Kirchhoff板越厚约束自动松弛以吸收横向剪切变形。2. 转角协调矩阵A_Delta_inv_DKMQ方法给出−3/2 · diag(1/(1φk))把边中点转角 Δβ 与角点转角 β 在高斯点上协调起来。3. 应变-位移矩阵 BB_b方法$$\mathbf{B}b \mathbf{B}{b(\beta)} \mathbf{B}{b(\Delta\beta)},\mathbf{A}\Delta^{-1},\mathbf{A}_u$$其中A_u由各边长与方向余弦dir_cos方法构成A_gamma、N_gamma负责转角到截面斜率的映射。4. 雅可比矩阵与高斯积分J方法用局部坐标构造 2×2 雅可比将参考坐标系求导转为物理坐标再在积分点上做BᵀDbB·det(J)加权求和得到 12×12 弯曲刚度子矩阵。5. 叠加膜力刚度弯曲之上再叠加一个等参数平面应力膜单元4 节点双线性 2×2 高斯积分扩张到 24×24 后直接相加即得单元总刚度ke ke_b ke_m。本构矩阵Dm面内与Db弯曲含t³/12因子在 Pynite/Plate3D.py 中实现还支持kx_mod/ky_mod正交各向异性刚度折减——这对模拟开裂混凝土非常实用。五、符号约定读懂内力结果的关键 有限元公式推导的最后一环是内力结果的提取与符号约定。PyNite 的梁弯曲符号约定如下y 向与 z 向各一张手绘图直观标注了荷载图、变形形状与 M(x)、V(x)、δ(x) 的关系板单元的结果提取同样讲究符号Pynite/Plate3D.py 的moment()方法通过 12 项系数矩阵C和曲率矩阵Q求得任意点(x, y)的Mx、My、Mxyshear()由弯矩对坐标求导合成Qx、Qymembrane()则在 4 个高斯点算应力后用外插形状函数H平移到目标位置。四边形板则直接在等参数坐标(ξ, η)上采样注意与矩形板的局部长度坐标区分开。六、实战检验与 Timoshenko 经典解对比 推导是否正确经典解说了算。示例 Examples/Rectangular Plate Bending - Qauds.py 用 1ft×1ft 的 Quad 网格建模一面 10ft×20ft 的四周固支墙体施加均布面压model.add_rectangle_mesh(MSH1, mesh_size, width, height, t, Concrete, 1, 1, [0, 0, 0], XY, element_typeQuad) model.analyze(check_staticsTrue)结果与 Timoshenko《板壳理论》表 35 的解析解对比弯矩幅值非常接近注意两者符号约定相反位移略偏大正是因为 DKMQ 考虑了横向剪切变形——厚板理论使然而非误差。PyNite四边形板单元剪壁结构有限元分析结果云图示例两个实用细节渲染的弯矩云图默认做平滑处理对汇聚于同节点的各单元角点应力取平均比直接用角点应力更准确网格对象还提供max_moment/min_moment方法直接提取极值。七、推导资料清单跟着官方一步步推 ✅想亲手复现整个公式推导按下面的路径走 Derivations/DMKQ Quad Element.ipynbDKMQ 四边形板单元的完整符号推导需 Jupyter Sympy运行全部单元格即可看到输出 Derivations/MITC4 Quad Element.ipynb对比学习 MITC4 另一种四边形板公式 Derivations/Rectangular Plate Element.ipynb12 项多项式矩形板的推导 Derivations/Fixed End Reactions - Linear Distributed Load.ipynb固定端反力推导理解等效节点力的基础 教材参考Katili (2015)、Bathe《Finite Element Procedures》、Logan《A First Course in the Finite Element Method》、Gallagher《Finite Element Analysis Fundamentals》完整书目见 Pynite/Quad3D.py 文件头注释。写在最后DKMQ 板单元的推导链条可以浓缩为一句话双线性插值位移 边中点二次转角 高斯点离散 Kirchhoff 约束 剪切耦合松弛 膜力叠加 弱弹簧稳定钻孔自由度。理解这条链条你不仅看懂了 PyNite 的板单元也拿到了分析一切等参数四边形板壳单元的万能钥匙 【免费下载链接】PyNiteA 3D structural engineering finite element library for Python.项目地址: https://gitcode.com/gh_mirrors/py/PyNite创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考
返回列表