二阶龙格库塔公式推导_有限元：推公式

在我本科的时候，学的东西都很浅很浅，但有时接触的东西很深，无奈，无法理解。

那会儿我们一帮理工科的人不用功，高数物理线代学的不明就里，曾以为只会好好做题便万事大吉，现在看，那时要是真的能做到好好做题也是不容易的。

那会儿我们一帮理工科的人不浪漫，搞了几年数模就错以为世间万物皆系于一模型之上，便幻想着用一组偏微分方程，来解出爱的模型。什么是爱呀？是建立爱的本构方程，找出女孩子的哈密顿量，解出它的本征值吗？

那会儿我们一帮理工科的人就像函数f(x)，g(x)，k(x)，每天都迷迷糊糊地输入着各种知识，输出却像死了一半，不见得有丝毫变化。

或许吧，曾经有个时间，我这个f(x)错以为找到了那个∫f'(x)dx，那个时候我知道，我唯一的答案是∫f'(x)dx。但那个时候的我还不知道，我只是∫f'(x)dx的一个选择。

当我步入研究生的时候，起初的那段日子相当难过，因为本科没有学过张量分析，还不明白△和▽是做什么用的，在学习有限元理论的路上，因为课本的公式推导没有过程，为此我也耗费了相当多的时间去手推。

言归正传，下面我将记录基本的一维度强形式和弱形式推导过程(不介绍FEM编程)。

好，

让我们定义一个u''(x)+l(x)=0

l(x) ∈[0,1]→R(x∈[0,1])

u''是二阶导数，du^2/dx^2

边界条件为u(1)=q，u'(0)=-h

这里u(1)=q一般叫做本质边界条件(第一类边界条件，essential boundary condition，Dirichlet boundary condition)

u'(0)=h叫做自然边界条件(第二类边界条件，natural BC，Newmann BC)

这就构成了一个强形式，一般来说是解不出来的，或者说不是很方便直接解出来。

实际上这里的u(x)为

其中y和z为哑指标

那么我们的下一步呢，则是得到弱形式。这里要提一下，推导弱形式是了解有限元方法基础中的基础，是每个了解有限元方法的同学都会的基础。

要推导弱形式，我们则要引入权重函数w，权重函数需要满足本质边界的齐次部分(homogeneous part)，即w(1)=0。

而构成弱形式的另一部分，即测试函数(trial solutions, candidate solutions)u仍然满足u(1)=q。

它还需要满足导数的平方可积，是因为FEM是一种Approximation方法，所以对于这种最小势能原理的微分方程，要保证变性能存在，而且强形式需要确保u二阶导的连续，推导弱形式添加了weight函数，由此降低了所需的导函数要求，即从要求二阶导连续，变成要求一阶导平方可积。

让我们回看强形式：

u''(x)+l(x)=0

两边乘以w

u''w+lw=0

再积分(积分上限1，下限0)，别忘了边界条件u'(0)=-h

∫u''wdx+∫lwdx=0

其中，∫u''wdx=wu'|-∫u'w'dx，由分部积分法得出。

整理一下便有了弱形式

∫u'w'dx=∫lwdx+w(0)h

到这里我们也看出来为什么定义第二类边界时用负号了，这里是为了方便构成弱形式。

弱形式的解满足强形式的解，强形式的解也满足弱形式的解，证明过程仍然是分部积分法，不再叙述。

至此，我们得到了弱形式，但这还不够，下一步是得到弱形式的解，有一种方法叫做伽辽金(Galerkin)方法。

简言之，假如空间被离散，特征长度为h

即我们的弱形式变成了：

其中

仍然满足边界条件：

由上面的式子，我们得到：

这边是Galerkin形式的方程

接下来求解过程

需要引入形函数(shape function)和定义：

NA便是形函数，NA(1)=0满足边界条件，即Wh(1)=0。

CA是已知的常数。

dA是常数，为满足边界条件，uh(1)=q，Nn+1(1)为1。

代入可得：

提出∑CA，发现∑CA可以被消去。

即

这样我们可以组建矩阵K，和F。

用矩阵来看

这样根据Kd=F，在matlab里便可求解出d=K-1F。

有了d之后，Galerkin方法的解便是：

~~~ 有时候会发一点数学建模， 和别的什么的。 可否，关注一下？