关于常系数高阶线性方程中所使用的辅助方程法

January 21, 2022

Abstract

这是一篇小报告，在上交之后我发表在这里。
对于常系数高阶线性方程的第二类待定系数法解，其的计算实在过于麻烦。于是[1]中提及一种简便方法来得到结果，那就是辅助方程法。但不管是网络上还是其他书籍里都少有相关方法的介绍，所以我在这里做一个简单的解释和总结。

原理

简单来说这是对欧拉恒等式的灵活运用，当在实数域时，一个 $a+bi$ 可以看作 $Rank=2$ ，但在复数域内则有 $Rank=1$ 。

辅助方程法的一个出发点是叠加原理，如下所示：

定理 1 (叠加原理). 设 $L[y]$ 为线性微分算子，若：

$y_1(x)$ 是 $L[y]=f_1(y)$ 的解；
$y_2(x)$ 是 $L[y]=f_2(y)$ 的解；

则有 $y_1(x)+y_2(x)$ 是 $L[y]=f_1(x)+f_2(x)$ 的解。

通过叠加原理我们可以组合不同的方程来得到一个组合解，与下面的独立性原理组合就能导出辅助方程方法。

定理 2 (独立性). 对于 $a+bi(a,b\in\mathbb R)$ ，其中 $a$ 与 $b$ 线性无关。

对于常微分方程这门课来说，我们利用复函数的目的是导出实函数解。考虑第二类其实是由第一类转换而得，显然我们可以通过第一类待定系数法来简化计算。

辅助方程法的步骤与构造方法

我们以 $y''+y=\frac12\cos x$ 为例。

例 1. 求解 $y''+y=\frac12\cos x$ 。

解. 对于 $\frac12\cos x$ ，我们可以对 $\frac12e^{ix}$ 展开得到。那么我们就设一方程为 $y''+y=\frac12e^{ix}$ 。

首先求得特征根 $\lambda_{1,2}=\pm i$ ，也即 $i$ 是一重特征根，由第一类方法即得解 $\phi_1(x)=-\frac i4xe^{ix}=\frac x4\sin x-\frac i4x\cos x$ 由独立性，显然 $y''+y=\frac12e^{ix}=\frac12\cos x+\frac i2\sin x$ 在 $\mathbb R$ 和 $\mathbb C\backslash\mathbb R\cup \{0\}$ 上有分别的解。这样就可以得出 $\frac x4\sin x$ 是 $y''+y=\frac12\cos x$ 的特解，而 $\frac x4\cos x$ 是 $y''+y=\frac12\sin x$ 的特解。 ◻

对于同样的例子，如果我们使用经典的第二类待定系数法来解，则会得到以下过程：

由一重特征根设特解为 $\phi(x)=x(A\cos x+B\sin x)$ ，对式子求导得 $\phi''(x)=(2B-Ax)\cos x-(A+Bx)\sin x$ 。由是得 $\phi''+\phi=2B\cos x+(Bx-A-Bx)\sin x=\frac12\cos x$ ，由 $\cos x$ 与 $\sin x$ 的正交性可得方程组 $\begin{align} 2B&=\frac12\\ (-A)&=0 \end{align}$ 由是解得 $B=\frac14$ ，即有同上的解。

上面的方程是一个简单的例子，我们看一个稍微复杂的例子来进行分析。通过这个例子我们将更加明了地认识到辅助方程法的优越性。

我们以 $y''-2y'+2y=xe^x\cos x$ 为例。

例 2. 求解 $y''-2y'+2y=xe^x\cos x$ 。

解. 对于 $xe^x\cos x$ ，很显然地我们找到 $xe^{(1+i)x}$ ，他通过Euler等式即可化为 $xe^x\cos x+ixe^x\sin x$ 。

由是我们解方程 $y''-2y'+2y=xe^{(1+i)x}$ 。

求解特征根得 $\lambda_{1,2}=1\pm i$ ，由是 $1+i$ 是一重特征根，即得解的形式为 $\phi(x)=x(ax+b)e^{(1+i)x}$ ，代入得解为 $\begin{aligned} \phi(x)&=x\left(-\frac i4x+\frac14\right)e^{(1+i)x}\\&=\frac x4e^x\left(\cos x+x\sin x+i(\sin x-x\cos x)\right) \end{aligned}$ 而 $xe^x\cos x$ 对应的解即是上面根的实部。所以原方程特解为 $\phi(x)=\frac x4e^x(\cos x+x\sin x)$ ◻

对于该方程我们用第二类待定系数法去解：

由于 $1+i$ 是一重特征根，则我们设解为 $\phi(x)=xe^x[(Ax+b)\cos x+(Cx+d)\sin x]$ ，则有 $\begin{multline} \phi'(x)=e^x ((-A x^2 - b x + C x^2 + 2 C x + d x + d)\sin x + \\(A x^2 + 2 A x + b x + b + C x^2 + d x)\cos x) \end{multline}$ $\begin{multline} \phi''(x)=-2 e^x ((A x^2 + 2 A x + b x + b - 2 C x - C - d)\sin x + \\(-2 A x - A - b - C x^2 - 2 C x - d x - d)\cos x) \end{multline}$ 由是有 $\begin{multline} \phi''(x)-2\phi'(x)+2\phi(x)=\\ 2 e^x ((-2 A x - b + C)\sin x + (A + 2 C x + d)\cos x)=xe^x\cos x \end{multline}$ 可列方程组 $\left\{ \begin{aligned} -4A&=0\\ 2(-b+C)&=0\\ 2(A+d)&=0\\ 4C&=1 \end{aligned} \right.$ 由是得解。

求解过程的复杂程度分析

我们可以用简单的迭代知识求得对一个式子朴素的求导次数。通过这种方法我们可以进一步考虑辅助方程法的优越性。

对于第一类换元法

我们记有 $\psi_{p,k}=\frac{\mathrm d{\psi_{p,k-1}}}{\mathrm dx},$

对于项 $\psi_{m,0}=Q_{m,0}(x)e^{\mu x}$ ，通过求导的乘法公式可以得出我们需要做 $m$ 次求导操作与 $m-1$ 次加法，得到一个新的式子 $\frac{\mathrm d{\psi_{m,0}}}{\mathrm dx}=\psi_{m,1}=Q_{m,1}(x)e^{\mu x}.$

我们设 $Q_{m,0}(x)=\sum_{k=0}^m{a_kx^k}$ ，则有 $\begin{align*} Q_{m,0}(x)&=\sum_{k=0}^m{a_kx^k}\\ Q_{m,1}(x)&=\sum_{k=0}^{m-1}{(a_k+\mu a_{k+1})x^k}+a_mx^m\\ Q_{m,2}(x)&=Q_{m,1}'+Q_{m,1}=Q_{m,0}''+2Q_{m,0}'+Q_{m,0}\\ \cdots\\ Q_{m,k}(x)&=\sum_{j=0}^k\binom{k}{j}Q_{m,0}^{(j)}\\ \cdots \end{align*}$ 对于 $Q_{m,0}^{(j)}$ ，其的每一项的系数为 $\frac{k!}{(k-j+1)!},k=j,\cdots,m$ 。统共要做 $j(m-j+1)$ 次乘法。对于系数又要做 $\min\{j,m-j\}$ 次乘法。则求 $Q_{m,k}(x)$ 要做 $\min\{j,m-j\}m(m-j+1)$ 次乘法，合并 $\sum_{j=1}^mj(m-j)=\frac16m(1+m)(m-1)$ 次同类项。最坏的情况下我们需要计算 $\sum_{l=0}^{m}b_lQ_{m,l},b_l\not=0$ 。则我们要做 $\sum_{j=0}^m\min\{j,m-j\}m(m-j+1)$ 次乘法，合并 $\frac16m^3(1+m)(m-1)$ 次同类项。

最后解方程组，若通过高斯消元法则需要 $O(m^3)$ 的时间复杂度消元。

如使 $m=3$ ，则要做 $33+36+27$ 次基本操作，假设每次操作有 $.1\%$ 可能性出错。则 $P_{WRONG}=1-(.999)^{96}=.092$ 。

对于第二类换元法

当使用第二类换元法时，由于要设两组 $Q_{m,0}$ ，则 $m$ 看作翻倍，但 $Q_{m,k}\equiv0,k>\frac m2$ ，如是我们可以修正得如下表。

	第一类	第二类
乘法	$\sum_{j=0}^m\min\{j,m-j\}m(m-j+1)$	$\sum_{j=0}^m\min\{j,2m-j\}2m(2m-j+1)$
合并同类项	$\frac16m^3(1+m)(m-1)$	$\frac13m^3(1+2m)(2m-1)$
解方程	$O(m^3)$	$O(m^3)$

取 $m=3$ ，则要做 $162+315+27$ 次基本操作，同样假设每次操作有 $.1\%$ 可能性出错。

则 $P_{WRONG}=1-(.999)^{505}=.396$ ，容易犯错。

同方法在其他范围内的应用

由 $e^x$ 的级数定义与Taylor展开很容易证明 $\begin{equation} \frac{\mathrm d{e^{\omega x}}}{\mathrm dx}=\omega e^{\omega x}, x\in\mathbb{R}, \omega\in\mathbb C \end{equation}$ 所以我们在对 $e^{\omega x}$ 的运算中可以把他当作是实函数去计算以简化操作。

概率论的特征函数

参考[2]，实际上这是对原密度函数做Fourier变换，也即 $\phi(t)=\int_{-\infty}^{\infty}e^{itx}p(x)dx$

例 3. 求 $Exp(\lambda)$ 的特征函数.

Proof. 考虑[expon]，显然地有 $\frac{\mathrm d{\frac{e^{\omega x}}{\omega}}}{\mathrm dx}=e^{\omega x}, x\in\mathbb{R}, \omega\in\mathbb C$ 。那么就有不定积分 $\int e^{\omega x}\mathrm dx=\frac{e^{\omega x}}{\omega}+C$ 。

而 $Exp(\lambda)$ 的PDF为 $p(x)=\left\{\begin{aligned} \lambda e^{-\lambda x}&,x\ge 0,\\ 0&,x<0, \end{aligned}\right.$ 由是有 $\begin{align*} \phi(t)&=\int_0^{+\infty}e^{itx}\lambda e^{-\lambda x}\mathrm dx\\ &=\lambda\int_0^{+\infty}e^{(it-\lambda)x}\mathrm dx\\ &=-\frac\lambda{it-\lambda} \end{align*}$ ¹ 附上教材解法： $\begin{align*} \phi(t)&=\int_0^{+\infty}e^{itx}\lambda e^{-\lambda x}\mathrm dx\\ &=\lambda\left[\int_0^{+\infty}\cos{(tx)}e^{-\lambda x}\mathrm dx+i\int_0^{+\infty}\sin{(tx)}e^{-\lambda x}\mathrm dx\right]\\ &=\lambda\left(\frac\lambda{\lambda^2+t^2}+i\frac t{\lambda^2+t^2}\right)=-\frac\lambda{it-\lambda} \end{align*}$ ◻

正态分布是一个更好的例子，我们通过变换形式而避免将函数展开。

例 4. 求 $N(0,1)$ 的特征函数.

Proof. $N(0,1)$ 的PDF为 $p(x)=\frac1{\sqrt{2\pi}}e^{-\frac{x^2}2}, -\infty<x<+\infty,$ 由是其特征函数为 $\begin{align*} \phi(t)&=\int_{-\infty}^{+\infty}e^{itx}\frac1{\sqrt{2\pi}}e^{-\frac{x^2}2}\mathrm dx\\ &=e^{-\frac{t^2}2}\int_{-\infty}^{+\infty}\frac1{\sqrt{2\pi}}e^{-\frac{(x-it)^2}2}\mathrm d(x-it)\\ &=e^{-\frac{t^2}2} \end{align*}$ ² 同样附上教材解法： $\begin{align*} \phi(t)&=\int_{-\infty}^{+\infty}e^{itx}\frac1{\sqrt{2\pi}}e^{-\frac{x^2}2}\mathrm dx\\ &=\frac{1}{\sqrt{2\pi}} \int_{-\infty}^{+\infty} \sum_{n=0}^\infty \frac{(itx)^n}{n!} e^{-\frac{x^2}{2}}\mathrm dx\\ &=\sum_{n=0}^\infty\frac{(it)^n}{n!}\left[\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{+\infty}x^ne^{-\frac{x^2}{2}}\mathrm dx\right] \end{align*}$ 方括号内是标准正态分布的 $n$ 阶矩 $E(X^n)$ ，当 $n$ 为奇数时由对称性 $E(X^n)=0$ ；当 $n$ 为偶数时，取 $n=2m$ ，则有 $E(X^n)=E(X^{2m})=(2m-1)!!=\frac{(2m)!}{2^m\cdot m!}$ 由是可得 $\phi(t)=\sum_{m=0}^\infty\frac{(it)^2m}{(2m)!}\frac{(2m)!}{2^m\cdot m!}=\sum_{m=0}^\infty\left(-\frac{t^2}{2}\right)^m\frac1{m!}=e^{-\frac{t^2}2}$ ◻

Laplace变换具有类似形式

Laplace变换和Fourier变换在形式上类似，我们可以通过类似方法计算。

例 5. 计算 $f(t)=e^{at}$ 的Laplace变换。

Proof. $\begin{align*} \int_0^{+\infty}e^{-st}e^{at}\mathrm dt&=\int_0^{+\infty}e^{(a-s)t}\mathrm dt\\ &=\frac{1}{a-s}\int_0^{+\infty}e^{(a-s)t}\mathrm d(a-s)t\\ &=\frac{1}{s-a} \end{align*}$ ³ ◻

袁荣 常微分方程; 高等教育出版社, 2012;

茆诗松; 程依明; 濮晓龙 概率论与数理统计教程; 高等教育出版社, 2011;

维基百科柯西积分定理 — 维基百科, 自由的百科全书 2021.

(https://math.stackexchange.com/users/1321/george-lowther), G.L. Characteristic Function of the Normal Distribution.

$|e^{itx}|\le1$ ，而 $e^{-\lambda x}\rightarrow0, x\rightarrow+\infty$ ，由柯西积分定理[3]↩︎
更严谨地，参考[4]。大意：由Lebesgue控制收敛定理可交换求导与求期望过程，由是得 $\mathbb E(e^{tX})$ 是解析的，而 $e^{-\frac{x^2}2}$ 同样是解析的。由解析延拓的唯一性，二函数若于实轴全相等则于任意值相等。↩︎
同样由柯西积分定理[3]↩︎