用matlab进行有限元法编程.doc

用matlab进行有限元法编程 Jack Chessa 2002年10月3日 1 引言 本文的目的在于给用户在使用matlab编写有限元代码方面遇到的一些问题作一个简短的概述和指引。估计读者对于有限元方法的理论知识已经有了基本的认识，因此，本文主要关注有限元方法的实施方法。本文将使用一个用mablab编写的用于分析静态线性弹性力学问题的有限元代码示例，以给大家说明如何进行有限元法的编程。 1.1 符号 为清楚起见，我们在本文中采用以下符号；粗体斜体字体V表示一个空间向量，其等同于在一个点上的位移或速度的空间维度，粗斜体字体d表示一个向量或矩阵的维数是未知数在离散系统即像STI内斯矩阵的系统矩阵，大写标表示节点编号，而小写标一般指沿笛卡尔单位向量的向量组成部分。所以，如果d是节点未知数的系统向量，uI是节点I的一个位移向量，那么 uIi是节点I在方向i上的位移向量即uI*ei。MATLAB的语法通常与数学符号混合使用，本文对此增加了增加了明确的解释。本文用打字机字体，字体指示MATLAB的语法被使用。 2 有关编写matlab程序 matlab编程语言在说明如何进行有限元法编程方便十分有用，因为其允许我们进行非常迅速的编写数值计算，并且其具有广阔的prede NED数学库。另外，矩阵（疏与密），矩阵和许多线性代数工具极影呗定义了，开发人员可以完全集中在没有定义这些数据结构的算法的实现上。广泛的数学和图形功能进一步使开发者免于开发这些功能本身或是找到适当的现有的库。一个简单的二维有限元程序，在matlab中只需要几百行代码，而在Fortran或C++程序中可能需要几千行。 尽管matlab编程语言在数学函数方面的功能非常完整，其依然有一些有限元的有利于发展为独立的功能的特定任务。这些已经被程序化，可参看前面提到的网站。 通常有贸易?FF这种易于开发。由于matlab是一种解释性的语言，每一行代码由matlab的命令行解释器解释并在运行时被顺序执行，因此，其运行时间可能比编译的编程语言，如Fortran或C + +，要长得多。应当指出的是内置的MATLAB功能已经编制并且效率极高，因此应尽可能多的使用。记住，是matlab的解说性质使得程序放缓，应不惜一切代价避免程序中的循环，尤其是嵌套的循环，因为这些可以使MATLAB程序的运行时间比可能需要的时间要长得多。通常的for循环可以通过使用matlab的量化处理消除。例如，下面的matlab代码设置矩阵A的行和列至零，并在对角线上放一个元素 for i=1:size(A,2) A(n,i)=0; end for i=1:size(A,1) A(i,n)=0; end A(n,n)=1; 这种代码一定不能使用因为下面的代码 A(:,n)=0; A(:,n)=0; A(n,n)=0; 在三行解释代码中做的事情是一样的，如同nr+nc+1解释行，其中A是一个nr×nc维矩阵。人们可以很容易地看到，当处理大型系统时（通常是用有限元代码的情况下），这可以迅速增加巨大的超载。for循环有时候是不可避免的，但令人惊讶的是，有几次是这种情况。建议人们当开发完一个matlab程序后，回去看看如何/是否可以消除任何循环。通过实践，这将成为matlab编程的第二特性。 3 一个典型有限元程序的组成部分 一个典型的有限元程序包含以下三个部分 预处理部分 处理部分 后处理部分 在预处理部分定义了问题陈述的数据和结构会被定。这些措施包括有限元离散，材料的性能，解决参数等等。在处理部分，有限元素对象即STI内斯矩阵，计算力矢量等被计算，边界条件被实施，系统被解决。在后处理部分，处理部分的结果将被分析。这里应力可能被计算并且数据进行可视化。本文中，我们将主要关注处理部分。许多前期和后期处理操作已经在matlab中进行编程并且在网上有参考文献；有兴趣的人可以直接参看MATLAB脚本文件或在MATLAB命令行中键入帮助“function name，以获得有关如何使用这些功能的进一步信息。 4 matlab中有限元数据结构 在这里，我们讨论的是在有限元法中使用，尤其在示例代码中实现的数据结构。人们可以想出多种储存有限元程序数据的方式，这有些随意，而我们尽力使用最灵活而且最有利于matlab的结构。这些数据结构的设计可能取决于所使用的编程语言，但通常不会与这里列出的结构明显不同。 4.1 节点坐标矩阵 由于我们进行有限元法编程，因此 我们需要一种代表元素离散域的方式。这样做，我们定义一组节点和用一定方式连接这些节点的元素。节点的坐标存储在节点坐标矩阵中。这是一个简单的节点坐标的矩阵（可以想象一下）。这个矩阵的维是nn×sdim，其中nn是节点的数量，sdim是空间维数。所以，如