微信公众号：[易木木响叮当]
关注可了解更多的[有限元数值仿真技巧]。问题或建议，请公众号留言;
如果你觉得木木同学对你有帮助，欢迎赞赏。

手撕梁单元UEL子程序

今天木木给大家带来的是梁单元的UEL子程序，作为杆系单元重要的组成成员，该部分是得好好讲讲了，而且要全面的讲，尽木木的知识能力所及讲的透彻，让刚入手UEL的小白们看的“起劲”！参考文档和源文件可在公众号后台回复：梁单元UEL即可。

今天分享的内容主要有：梁单元介绍、等效节点荷载处理、状态变量的自定义输出、理论解对比等等，精彩超多，请进一步往下翻阅

模型讲解

梁单元的讲解在前几期的有限元基础编程的推文中有详细讲述，感兴趣的小伙伴可以跳转相关页面查阅。

案例简介

我们今天讲的案例如图二所示，梁长，分布荷载为-10，将左端点的U1、U2、UR3自由度束缚，右端点的U2自由度束缚。

理论解

这里我们简单讲述一下变形协调条件，图3 所示（图2分解）为一承受均布荷载的悬臂梁，在B端点产生的位移在数值上可等价于由支反力产生的位移，有：

由以上公式可得到右端点节点力、弯矩以及左端点节点力。

子程序讲解

编制矩阵操作函数

为了程序的易读性，我们可以将经常使用的的计算部分编写成一个子程序，使用时直接call一下即可，比如在UMAT中，我们会常常用到Mises应力公式，所以我们可以将Mises的求解编写为一个子程序，在UMAT中调用即可，大大减少了代码量，更能提升易读性。

本程序中，我编写了矩阵置零函数，用于数组的初始化，大家在自己的程序中编制矩阵相乘函数和矩阵求逆函数，大大提升编写代码的效率。

等效荷载处理

在本次案例中，我们注意到，梁体受到向下的均布荷载，以往的UEL案例分析中，我们没考虑过这种面力、体力的作用，可以参考一下势能泛函：

具体含义可见附带文档。本案例的面力处理方式：将其转化为单元节点荷载，转化的原理可参考附带文档，也可查阅有限元相关教材，推荐曾攀老师的《有限元基础教程》。

状态变量输出

为了与理论解析解做对比，可以以状态变量的形式将右端点的节点力和弯矩以及左端点的节点力输出到指定位置，这里输出到指定文件中。涉及到Fortran语言中文件的使用，编制了计算节点力的子程序，如下所示。

状态变量输出

即可在'D:\TEMP'位置处创建STRESS.DAT文件，按照以上语句进行写出。

理论解析解对比

与解析解一致。最后来看一下位移云图吧！

今天的分享就止于此了，下期我们再见啦~~喜欢本次案例的小伙伴可以点点小赞和在看哦，欢迎转发至朋友圈，让更多的对二次开发感兴趣的伙伴看到。