常微分方程解法总结(实用13篇)

常微分方程解法总结 第1篇

到目前为止，我们讨论了步长取常数的常微分方程数值解法。

但实际中，方程的曲线在有的地方变化缓慢，有的地方变化剧烈。在变化缓慢的地方可以用大一点的步长，在变化剧烈的地方可以取小一点的步长。于是就有了自动调整步长的算法。

自动调整步长的算法可以避免保证计算精度的同时提高效率。

自适应步长控制((adaptive step-sixe control)的方法在执行过程中需要估计每一步的局部截断误差。然后，根据该误差估计值确定步长应该变大还是缩小。

在单步法中加入自适应步长控制主要有两种方式： 第一种方式是，计算同阶龙格-库塔法在两个不同步长F的预测值，将这两个值之差作为误差的估计值； 第二种方式是，算出两个不同阶龙格-库塔法的预测值，将它们的差作为局部截断误差的估计值。

步长控制的第一种方式，就是先用步长h计算一个 y h y_h yh​，然后再用步长为h/2分两步计算 y h y_h yh​，比较这两个计算结果的差值，将这个差值作为局部截断误差的估计值。 若记 y 1 y_1 y1​为单步预测值， y 2 y_2 y2​为两个半步所得到的预侧值，则误差 Δ \Delta Δ可表示为: Δ = y 2 − y 1 \Delta = y_2-y_1 Δ=y2​−y1​（）

这个 Δ \Delta Δ就可用于作为控制步长的标准。

此外， Δ \Delta Δ还可以用于校正 y 2 y_2 y2​的预测值。譬如说，对于四阶龙格-库塔法，其估值原本具有四阶精度。如果用式()校正后，可达到五阶精度。// 具体原因，可以结合“龙贝格积分的递推思路”？ y 2 ← y 2 + Δ / 15 y_2 \gets y_2 + \Delta / 15 y2​←y2​+Δ/15（）

步长控制的第二种方式，先用四阶龙格-库塔法计算一个预测值 y 1 y_1 y1​，再用五阶龙格-库塔法计算一个预测值 y 2 y_2 y2​，（当然，也可以考虑用一个二阶+一个三阶；也或者用一个二阶+一个四阶，都是可以的。）然后比较这两个计算结果的差值。

该方法的一个缺点是，计算开销增长过大。例如，若使用四阶和五阶格式作为预测步，则每步总共要计算10次函数值（四阶4次，五阶6次）。

有没有什么方法可以客服上述问题呢？？

于是我们就思考，如何重复利用计算过的函数值。譬如说，如果在五阶龙格-库塔法中能使用所有四阶龙格-库塔法所计算的函数值，这样只需要计算6次函数值，就得到一个误差估计。

沿着这个思路，因为四阶和五阶的龙格-库塔法都有很多种不同的系数组合方式，如果使用下列的四阶和五阶公式： y i + 1 = y i + ( 37 378 k 1 + 250 621 k 3 + 125 594 k 4 + 512 1771 k 6 ) ∗ h y_{i+1}=y_i + (\frac{37}{378}k_1 + \frac{250}{621}k_3 + \frac{125}{594}k_4 + \frac{512}{1771}k_6)*h yi+1​=yi​+(37837​k1​+621250​k3​+594125​k4​+1771512​k6​)∗h（） y i + 1 = y i + ( 2825 27648 k 1 + 18575 48384 k 3 + 13525 55296 k 4 + 277 14336 k 5 + 1 4 k 6 ) ∗ h y_{i+1}=y_i + (\frac{2825}{27648}k_1 + \frac{18575}{48384}k_3 + \frac{13525}{55296}k_4 + \frac{277}{14336}k_5 + \frac{1}{4}k_6)*h yi+1​=yi​+(276482825​k1​+4838418575​k3​+5529613525​k4​+14336277​k5​+41​k6​)∗h（） 其中， 利用式()和式()计算2个预测值，然后用两个预测值之差作为误差的估计值。

思考，在自适应龙格-库塔法或步长——对分法，因为一阶和二阶是从一个点到两个点来计算预测值，然后用两个预测值只差计算出误差估计值。然后误差估计值又可用于改进预测值，能提高一阶精度。

思路一：是不是也可以用二阶和四阶的两个预测值来计算估计值，然后用该估计值来改进预测值，会不会带来更好的效果？？

思路二：因为之前提到过，五阶是不太划算的，因为五阶要计算六次函数值。那么是否可以考虑用三阶和四阶搭配来计算两次预测值？或者二阶和四阶搭配？能否在各种搭配中找出更高效的方法。

上一节介绍了估计局部截断误差的方法，利用它可以调整步长。若误差很小，则增大步长；若误差很大，则减小步长。为此，Press等人提出了下面的标准: h n e w = h p r e s e n t ∣ Δ n e w Δ p r e s e n t ∣ α h_{new} = h_{present}|\frac{\Delta_{new}}{\Delta_{present}}|^{\alpha} hnew​=hpresent​∣Δpresent​Δnew​​∣α（） 为啥 α \alpha α取或？因为四阶和五阶开4次方和开5次方正好是1/4和1/5？

式()中的关键参数是 Δ n e w \Delta_{new} Δnew​，因为它可用于指定期望精度。

一种方法是，将 Δ n e w \Delta_{new} Δnew​与一个相对误差水平联系起来。但是当解越过零点时，会产生一些问题。处理这种情况的一个更通用的方法是，将 Δ n e w \Delta_{new} Δnew​定义为 Δ n e w = ε y s c a l e \Delta_{new} = \varepsilon y_{scale} Δnew​=εyscale​ () 其中， ε \varepsilon ε为整体容许水平。 y s c a l e y_{scale} yscale​的选取将决定误差的尺度。其中常见 y s c a l e y_{scale} yscale​采用下式： y s c a l e = ∣ y ∣ + ∣ h d y d x ∣ y_{scale} = |y|+|h\frac{d_y}{d_x}| yscale​=∣y∣+∣hdx​dy​​∣

图中的点表示由自适应步长例程得出的预测值

常微分方程解法总结 第2篇

则方程有特解 y ∗ = x Q m ( x ) y^{*}=xQ_{m}(x) y∗=xQm​(x)

常微分方程解法总结 第3篇

假设这对解是 α ± i β \alpha\pm i\beta α±iβ

此时方程的通解是

Y ( x ) = e α x ( C 1 c o s ( β x ) + C 2 s i n ( β x ) ) Y(x)=e^{\alpha x}(C_{1}cos(\beta x)+C_{2}sin(\beta x)) Y(x)=eαx(C1​cos(βx)+C2​sin(βx))

对于 f ( x ) f(x) f(x)和特征根的情况，对特解的情况做如下归纳：

常微分方程解法总结 第4篇

则方程有特解 y ∗ = Q m ( x ) y^{*}=Q_{m}(x) y∗=Qm​(x)

常微分方程解法总结 第5篇

举例：假设跳伞人的下落速度v于时间有如下关系：

d v d t = g − c m v \frac{dv}{dt} = g-\frac{c}{m}v dtdv​=g−mc​v ()

其中g为重力常数，m为质量，c为阻力系数。

被微分的量v因变量，与v有关的变量t称为自变量。 如果函数只有一个自变量，那么方程就称为常微分方程（ordinary differential equation，ODE）。 如果函数含有两个或者多个自变量，则成为偏微分方程（partialdifferential equation，PDE）。

此外，微分方程也可以根据阶数来分类：最高阶导数是一阶导数，则称为一阶微分方程（first-order-equation）；最高阶导数是二阶导数，则称为二阶微分方程（second-order-equation）。例如如（）中就是一阶微分方程，下式（）就是一个二阶微分方程。

m d 2 x d t 2 + c d x d t + k x = 0 m\frac{d^2x}{dt^2} +c\frac{dx}{dt} + kx = 0 mdt2d2x​+cdtdx​+kx=0 ()

高阶微分方程能简化成一阶方程组。考虑上式（），定义新变量y，令 y = d x d t y=\frac{dx}{dt} y=dtdx​ ()

对上式取微分得： d y d t = d 2 x d t 2 \frac{dy}{dt}=\frac{d^2x}{dt^2} dtdy​=dt2d2x​ ()

将式（）和（）代入式（）中得到： m d y d t + c y + k x = 0 m\frac{dy}{dt} +cy + kx = 0 mdtdy​+cy+kx=0 ()

于是原来的二阶微分方程（）可等价于两个一阶方程组（）和（）。

同样的，其他的n阶微分方程也可以用类似的方式简化。

常微分方程通常采用解析积分得方法来求解。如对于式()，先乘以dt，在进行积分得到： v = ∫ ( g − c m v ) d t v = \int{(g-\frac{c}{m}v)dt} v=∫(g−mc​v)dt ()

对于上式（），是可以精确的推导出该积分得函数表达式的。因为该方程是线性的。

但在实际中，很多方程（是非线性的）精确解是无法求出的。于是提出了一个方法，就是将方程线性化。 （ n n n阶）线性常微分方程的一般形式是： a n ( x ) y ( n ) + . . . + a 2 ( x ) y ( 2 ) + a 1 ( x ) y ′ + + a 0 ( x ) y = f ( x ) a_n(x)y^{(n)}+...+a_2(x)y^{(2)}+a_1(x)y'++a_0(x)y = f(x) an​(x)y(n)+...+a2​(x)y(2)+a1​(x)y′++a0​(x)y=f(x) () 其中， y ( n ) y^{(n)} y(n)是y关于x的n阶导数， a n ( x ) a_n(x) an​(x)和 f ( x ) f(x) f(x)都是关于x的函数。因为该方程中未出现因变量y与其导数的乘积，也没有出现非线性函数。所以认为它是线性的。

如下式（）是一个非线性微分方程： d 2 x d t 2 + g l s i n ( x ) = 0 \frac{d^2x}{dt^2} +\frac{g}{l} sin(x)= 0 dt2d2x​+lg​sin(x)=0 () 由于含有 s i n ( x ) sin(x) sin(x)为非线性函数，故该微分方程是非线性的。

线性常微分方程是可以通过解析法求解的。但是，大部分非线性方程无法精确求解。

拿一个简单的方程举例。首先给定函数： y = − x 4 + 4 x 3 − 10 x 2 + x + 1 y= y=−x4+4x3−10x2+x+1 （） 这是一个四次多项式。对其进行微分，就得到一个常微分方程： d y d x = − 2 x 3 + 12 x 2 − 20 x + \frac{dy}{dx}=-2x^3+12x^2-20x+ dxdy​=−2x3+12x2−20x+ （）

对式（）乘以dx，在进行积分得到： y = ∫ ( − 2 x 3 + 12 x 2 − 20 x + ) d x y=\int{(-2x^3+12x^2-20x+)}dx y=∫(−2x3+12x2−20x+)dx （） 应用积分法则得出解为： y = − x 4 + 4 x 3 − 10 x 2 + x + C y= y=−x4+4x3−10x2+x+C （） 除了相差一个C外，其余都与原函数相同。这个C称为积分常数（constant of integration）。 出现一个任意常数C表明，积分的结果并不算是唯一的。无限多个常数C对应无限多个可能的函数，都满足微分方程。下图给出了6个满足条件的函数：

为了将解完全确定下来，微分方程通常伴随有辅助条件（auxiliary conditions）。对于一阶常微分方程，有一类被称为**初值（initial value）**的辅助条件，这类条件用于确定常数值，从而使得解是唯一的。例如，给式（）添加初始条件x=0，y=1。带入式（）中，可推导出C=1。于是就得到了唯一解： y = − x 4 + 4 x 3 − 10 x 2 + x + 1 y= y=−x4+4x3−10x2+x+1 这个解同时满足常微分方程和指定的初始条件。

当处理n阶微分方程时，就需要n个条件来确定唯一解。如果所有的条件都是在自变量同一值处指定的，那么问题就称为初值问题（initial-value problem）。与之相对的，边值问题（boundary-value problems），就是指在自变量的不同值处指定初始条件。

常微分方程解法总结 第6篇

T(p)=\int_0^{+\infty}f(t)e^{-pt}dt

f(t)=\frac1{2\pi j} \int_{\beta-j\infty}^{\beta+j\infty}T(p)e^{pt}dp

感觉有点像fly？

有没有觉得好像欧拉方程（完整版）（bushi

用来解决高阶线性很好玩，先做个变换：

s^2X(s)-2sX(s)+2X(s)= \mathcal{L} [2e^tcost]=\frac{2(s-1)}{(s-1)^2+1}

X(s)=\frac{2(s-1)}{[(s-1)^2+1]^2}

x(t)=te^tsint

x^2y''+xy'+(x^2-v^2)y=0

要用广义幂级数，感觉好麻烦

不是未完待更新，而是早就考完了

常微分方程解法总结 第7篇

则方程有特解 y ∗ = e α x ( A 1 c o s ( β x ) + A 2 s i n ( β x ) ) y^{*}=e^{\alpha x}(A_{1}cos(\beta x)+A_{2}sin(\beta x)) y∗=eαx(A1​cos(βx)+A2​sin(βx))

常微分方程解法总结 第8篇

则方程有特解 y ∗ = x k Q m ( x ) y^{*}=x^{k}Q_{m}(x) y∗=xkQm​(x)

其中， Q m ( x ) = b 0 + b 1 x + … … + b m x m Q_{m}(x)=b_{0}+b_{1}x+……+b_{m}x^{m} Qm​(x)=b0​+b1​x+……+bm​xm，

b i ( i = 0 , 1 , … … m ) b_{i}(i = 0,1,……m) bi​(i=0,1,……m)是需要带回原方程来确定的系数。

常微分方程解法总结 第9篇

设 λ \lambda λ是特征方程的一个1重实数解。

则应该为它在通解中增加这样一项：

C e λ x Ce^{\lambda x} Ceλx

常微分方程解法总结 第10篇

设 α ± i β \alpha\pm i\beta α±iβ是特征方程的一对1重共轭复数解。

则应该为它在通解中增加这样一项（如果拆开括号就是2项）：

e α x ( C 1 c o s ( β x ) + C 2 s i n ( β x ) ) e^{\alpha x}(C_{1}cos(\beta x)+C_{2}sin(\beta x)) eαx(C1​cos(βx)+C2​sin(βx))

常微分方程解法总结 第11篇

设 α ± i β \alpha\pm i\beta α±iβ是特征方程的一对k重共轭复数解。

则应该为它在通解中增加这样一项（如果拆开括号就是2k项）：

e α x [ ( C 1 + C 2 x + . . . + C k x k − 1 ) c o s ( β x ) + e^{\alpha x}[(C_{1}+C_{2}x+...+C_{k}x^{k-1})cos(\beta x)+ eαx[(C1​+C2​x+...+Ck​xk−1)cos(βx)+

( D 1 + D 2 x + . . . + D k x k − 1 ) s i n ( β x ) ] (D_{1}+D_{2}x+...+D_{k}x^{k-1})sin(\beta x)] (D1​+D2​x+...+Dk​xk−1)sin(βx)]

当你为每一个特征方程的解写出对应的项后，把他们加起来，就得到齐次方程的通解了。

常微分方程解法总结 第12篇

假设该解等于 r r r，

此时方程的通解是

Y ( x ) = ( C 1 + C 2 x ) e r x Y(x)=(C_{1}+C_{2}x)e^{rx} Y(x)=(C1​+C2​x)erx.

常微分方程解法总结 第13篇

x^{(n)}+a_1(t)x^{(n-1)}+a_2(t)x^{(n-2)}+\cdots+a_n(t)x=f(t)

向量形式：

\frac{d\overrightarrow{x}}{dt}=A\overrightarrow{x}+\overrightarrow{F}(t)

\begin{pmatrix} 0 & 1 &\cdots &0\\ 0 & 0 &\cdots &0\\ \vdots & \vdots &\ddots &\vdots \\ -a_n(t) & -a_{n-1}(t) &\cdots &-a_1(t)\end{pmatrix}

若函数组xi线性相关，则W(t)=0。(但是W(t)！=0不一定函数族线性无关)

高阶线性微分方程线性无关的解组的Wronsky行列式永不为0.（如果行列式有零点那就线性相关）

W(t)=W(t_0)exp(-\int_{t_0}^ta_1(s)ds)

自己看上面的方程，a1是次高导数项的系数

用途：一个k阶方程，已知了k-1个解，还差一个解，就把W(t)列出来，右边等于一个常数C乘上exp积分。

y^{(n)}+a_1y^{(n-1)}+\cdots+a_ny=0

跟之前一阶一样的，把导数换成次方解方程,拿到特征根塞给exp。如果有重根，就在exp前面堆多项式，如果有复根，就写成\rho e^{i\theta}然后分开设expsin+expcos，和组合数学那里一样的。

x^{(n)}+a_1(t)x^{(n-1)}+a_2(t)x^{(n-2)}+\cdots+a_n(t)x=f(t)

ODE ppt讲这里没有一图流，直接看组合数学得了

常数部分=r^nb(n)

特解形式=n^m(k_0+k_2n+\cdots+k_qn^q)r^n

可以拆成三个部分相乘，第一部分：n的{重根}次方，第二部分：和等式右侧多项式次数相等的未定系数多项式，第三部分：直接把特征根部分搬下来

先把几个通解基搞出来，再把前面系数写成c(t)带进去解c（t）.等式常数部分怪怪的用不了凑特解就用这个。

设线性通解为y=c_1 x_1+\cdots+c_n x_n，则常数变易法得到的方程组为；

y''-y=xe^xcosx

y''+y=4sinx

用常数变易

y''+y=\frac{1}{sin^3x}

就不能叫什么欧拉第一方程 第二方程？还有个变分也叫欧拉方程

x^{n}y^{(n)}+a_1x^{n-1}y^{(n-1)}+\cdots+a_{n-1}xy'+a_ny=f(x)

特点：x的次数和y的导数次数一致

更新解法：更适合物理系思路的理解方式

令x=e^t

\dot y=y'\dot x=y'x

\ddot y=(y'x)'\dot x=(y''x+y')x=y''x^2+y'x

...

D(D-1)...(D-(n-1))y=x^ny^{(n)}

下面是原ppt讲法

令D=\frac{d}{dt}，则：

xy'=Dy

x^2y''=D(D-1)y

...

代回原方程，化为g(D)*y=\varphi(x)的形式，相当于是个非齐次线性方程，用特征方程随便解下就行了。记得最后把t换成x。