【求助】关于一维热传导方程的有限元解
毕业论文中要用到，但是没学过有限元。这里达人挺多，来试试！题目见附件。 这个问题存在精确解，为什么还有求有限元解？ 嘿嘿，为了说明问题的，并不是想说它有多难~~ Discrete domain: 1--------------------------2----------------------------3
x1 Element 1 x2 Element 2 x3
Assume we are using linear element and the element length is h, the shape functions are
For element 1, N1(x) = (x2-x1)/h, N2(x) =(x-x1)/h
For element 2, N2(x) = (x3-x)/h, N3(x) = (x-x2)/h
Take an arbitrary element xi--------------------------xj and apply the variational form, we have:
int du/dx dv/dx = int 2 * v(x)
Approximate u by u = ui*Ni(x) + uj*Nj(x), and let v = Ni(x) and Nj(x) respectively, we have
int <dNi(x)/dx*dNi(x)/dx, dNj(x)/dx*dNi(x)/dx> * <ui, uj>^T = int 2*Ni(x) for v=Ni(x)
int <dNi(x)/dx*dNj(x)/dx, dNj(x)/dx*dNj(x)/dx> * <ui, uj>^T = int 2*Nj(x) for v=Nj(x)
Therefore you can plug the parameters from Element 1 and Element 2 into the above formula, you will get the element stiffness matrix and load vector.
Now we need to apply the boundary conditions:
For u(0) = 0, we can put something like u(0) - epsilon*du/dx(0) = 0, epsilon is a very small number,
If we apply the variational form for the element touching the left boundary, we have
int du/dx * dv/dx = du/dx(h)*v(h) - du/dx(0)*v(0)
We know du/dx(0) = u(0)/epsilon, then we have:
int du/dx * dv/dx = du/dx(h)*v(h) - 1/epsilon*u(0)*v(0). We can forget about du/dx(h)*v(h)since it will be cancelled with the contribution from another element when assembling the element stiffness matrix. So we have actually:
int du/dx * dv/dx + 1/epsilon*u(0)*v(0) = 0
Note v(0) = N1(0) = 1, v(0) = N2(0) = 0. But noting u(0) = N1(0)*u1 + N2(0)*u2, we can get the element stiffness matrix by including the boundary conditions. This is the penalty methods. Same way, we can enforce u(1)=1. 谢谢tonnyw
的热心帮助，十分感谢～～
页:
[1]