规划有限元法

Question

我试着教自己如何使用有限元法。

我的所有代码都是从以下链接页16-20 http://homepages.cae.wisc.edu/~suresh/ME964Website/M964Notes/Notes/introfem.pdf改编而来的。

我正在Matlab中编程，对单个8节点立方体单元进行有限元分析。我已经定义了xi，eta，zeta局部轴(我们现在可以把它想成x，y，z)，所以我得到了以下形状函数：

%%shape functions

zeta = 0:.01:1;
eta = 0:.01:1;
xi = 0:.01:1;

N1 = 1/8*(1-xi).*(1-eta).*(1-zeta);
N2 = 1/8*(1+xi).*(1-eta).*(1-zeta);
N3 = 1/8*(1+xi).*(1+eta).*(1-zeta);
N4 = 1/8*(1-xi).*(1+eta).*(1-zeta);
N5 = 1/8*(1-xi).*(1-eta).*(1+zeta);
N6 = 1/8*(1+xi).*(1-eta).*(1+zeta);
N7 = 1/8*(1+xi).*(1+eta).*(1+zeta);
N8 = 1/8*(1-xi).*(1+eta).*(1+zeta);

[N]矩阵按我正在读的课文排列如下：

%N Matrix
N= [N1 0 0 N2 0 0 N3 0 0 N4 0 0 N5 0 0 N6 0 0 N7 0 0 N8 0 0;
    0 N1 0 0 N2 0 0 N3 0 0 N4 0 0 N5 0 0 N6 0 0 N7 0 0 N8 0;
    0 0 N1 0 0 N2 0 0 N3 0 0 N4 0 0 N5 0 0 N6 0 0 N7 0 0 N8];

要找到[B]矩阵，我必须使用以下[D]矩阵：

%%Del Matrix for node i
%[  d/dx   0     0
%    0    d/dy   0
%    0     0    d/dz        . . .
%   d/dy  d/dx   0
%    0    d/dz  d/dy
%   d/dz   0    d/dx   ]

它是一个在[N]上运行的操作员。(B=DN)

稍后，如本文所示，我将计算这个[B]矩阵在这个元素的体积上的积分。

因此，我的问题是，如何将这些多项式形状函数存储在一个矩阵中，对它们进行微分运算，然后对它们进行数值积分。我可以用我现在设置的方式来判断，它不能工作，因为我已经将函数定义为区间[0,1]上的向量，然后将这些向量存储在[N]矩阵中。然后利用diff()函数进行适当的判别，找出[B]矩阵。但是，由于[B]的矩阵元现在是区间[0,1]上的向量，我认为这会引起问题。你们怎么做这些计算，在我上面贴出的教科书里描述的？

Answer 1

用匿名函数解决了我的问题，并将多项式存储在一个符号矩阵中。示例：

syms xi eta zeta
N1= ... %type out in terms of xi eta and zeta
.
.
.
dN1dXi = diff(N1,xi) %symbolic differentiation with respect to xi

还可以在需要时执行符号集成：

intN1 = int(N1,xi,lowerLimit,upperLimit) %symbolic integration with respect to xi

当准备在实际值中替换以评估符号函数时：

subs(N1,{xi,eta,zeta},{value1,value2,value3})

Answer 2

关于如何从参数域(0,1^)映射到物理域，您应该检查第24页。

Answer 3

虽然我觉得你可以按你说的做，用象征性的。我认为Matlab中的符号计算非常耗时。

我会手动选择派生N，并将其存储为dN，并在需要时使用它。

致以敬意，

德语