基本解和 Green 函数的区别
如果将微分方程写成抽象的算子形式: \[
L[u] = f \text{ on } \Omega
\]
特别的, 当函数 \(f\) 为 Dirac 函数
\(\delta\) 时, 微分方程的解称为算子
\(L\) 的基本解: \[
L [G] = \delta
\]
注意到上述过程中, 并没有提到方程需要满足什么边界条件, 这也是基本解和
Green 函数之间的区别. Green 函数需要考虑边界的影响,
基本解只需要考虑函数在区域内满足什么微分方程. 在下文中,
我还是按照 Mura 书中的表述, 将以求解基本解方式得到的函数称为 Green 函数.
当微分方程右端项为一般的函数 \(f\) 时,
此时微分方程的解 \(u\) 可以表示为Green
\(G\) 与 \(f\) 的卷积:
\[
u = G*f = \int_{-\infty}^{\infty} G(x-x^{\prime}) f(x^{\prime})
\mathrm{d} x^{\prime}
\]
张量值函数的卷积定义
对于张量值方程 \[
\boldsymbol{L}[\boldsymbol{u}] = \boldsymbol{f}
\] 应用 Fourier 变换, 方程组转化为关于 \(\hat u\) 的代数方程: \[
\hat L_{ij}(\xi) \hat u_{j} = \hat f_i
\] 对方程系数矩阵 \(\hat
L_{ij}\) 求逆, 再作 Fourier 逆变换, 就得到解 \(u\) 表示为 \[
u_i = \mathscr{F}^{-1} [\hat{L}_{ij}^{-1} \hat{f}_j]
\] 定义张量值函数之间的卷积为 \[
A * f: L(\Omega,\mathbb{R}^{n\times n}) \times L(\Omega,\mathbb{R}^{
n}) \mapsto L(\Omega,\mathbb{R}^{n}) \\
A * f \triangleq A_{ij} * f_{j}
\] 上式默认了指标求和约定. 这也就是说,
张量的代数运算保存在卷积的积分号内.
张量值方程的基本解定义为 \[
G_{ij} = \mathscr{F}^{-1} [\hat{L}_{ij}^{-1}]
\] 因此, 方程解可以表示为 \[
\boldsymbol{u} = \boldsymbol{G} * \boldsymbol{f}
\]
含本征应变弹性方程的基本解
对于含本征应变的控制方程: \[
L_{ijkl}~u_{k,lj} - L_{ijkl}~\mu_{kl,j} = 0
\]
我们希望能够找到该方程的 Green 函数 \(G\), 使得方程的解可以表示为 Green
函数与本征应变的卷积. 这需要借助于 Fourier 分析方法. 对方程进行 Fourier
变换, 得到
\[
L_{ijkl}\xi_{l}\xi_{j}~\hat{u}_{k} + i \xi_{j}L_{ijkl}~\hat{\mu}_{kl} =
0
\]
记 \(K_{ik}(\xi) \triangleq
L_{ijkl}\xi_{l}\xi_{j}\), \(X_i
\triangleq -i \xi_{j}L_{ijkl}~\hat{\mu}_{kl}\), 就得到关于 \(\hat{u}_k\) 的线性方程:
\[
K_{ik} \hat{u}_k = X_i
\]
如果 \(L_{ijkl}\) 具有主对称性,
也就是说, \(L_{(ij)[kl]} =
L_{[kl](ij)}\), 那么矩阵 \(K_{ik}(\xi)\) 是对称矩阵. 再应用 Fourier
逆变换, 就可以得到解 \(u(x)\):
\[
u_{i}(x) = \mathscr{F}^{-1} [\hat{u}_{i}(\xi)]
= \mathscr{F}^{-1}[ -i\xi_{j}~ L_{kjmn}~K_{ik}^{-1}(\xi)\hat{\mu}_{mn}]
\]
将三维 Fourier 逆变换的公式代入得到
\[
\begin{aligned}
u_{i}(x)
&= - \frac{1}{(2\pi)^3}
\int i\xi_{j} L_{kjmn} K_{ik}^{-1}(\xi) e^{i\xi \cdot x} \mathrm{d} \xi
\int \mu_{mn}(x^{\prime}) e^{-i\xi \cdot x^{\prime}} \mathrm{d}
x^{\prime}\\
&= - \int L_{kjmn}
\frac{\partial}{\partial x_j} \left( \frac{1}{(2\pi)^3}
\int K_{ik}^{-1}(\xi) e^{i\xi \cdot (x-x^{\prime})}
\mathrm{d} \xi \right)
\mu_{mn}(x^{\prime}) \mathrm{d} x^{\prime}
\end{aligned}
\]
定义 Green 函数为弹性方程逆矩阵的 Fourier 逆变换:
\[
G_{ik}(x) \triangleq \mathscr{F}^{-1}[K_{ik}^{-1}](x)
=\frac{1}{(2\pi)^3} \int K_{ik}^{-1}(\xi) e^{i\xi \cdot x} \mathrm{d}
\xi
\]
那么位移 \(u_i\) 可以表示为 Green
函数与本征应变的卷积
\[
u_{i}(x) = -\int L_{jkmn} G_{ik,j}(x-x^{\prime}) \mu_{mn}(x^{\prime})
\mathrm{d} x^{\prime}
\]
式中, Green 函数下标 \(\square_{,j}\) 表示对变量 \(x_j\) 求偏导 (注意区别变量 \(x\) 与 \(x^{\prime}\)). 或者将上式记为卷积的形式:
\[
u_i = G_{ik,j} * \left( -L_{jkmn} \mu_{mn} \right)
\] ## 验证 Green 函数所满足的微分方程
对 Green 函数求二阶导数, 得到 \[
G_{ij,kl}(x) = \frac{1}{(2\pi)^3}
\int -K_{ij}^{-1}(\xi) \xi_k \xi_l e^{i\xi \cdot x} \mathrm{d} \xi
\]
方程两端乘 \(L_{jknl}\),
等式依旧成立, 得到
\[
\begin{align}
L_{jknl}G_{ij,kl}(x)
&= -\frac{1}{(2\pi)^3}
\int K_{ij}^{-1}(\xi) L_{jknl} \xi_k \xi_l e^{i\xi \cdot (x-x^{\prime})}
\mathrm{d} \xi
\end{align}
\]
注意到 \(K_{jn}(\xi) :=
L_{jknl}\xi_{k}\xi_{l}\), 所以有 \[
L_{jknl}G_{ij,kl}(x)
= -\delta_{in} \frac{1}{(2\pi)^3}\int e^{i\xi \cdot x} \mathrm{d} \xi
= -\delta_{in} \delta(x)
\] 整理一下指标, 得到 \[
L_{ijkl}G_{nk,lj}(x) + \delta_{in} \delta(x) = 0
\]
展开指标 \(i\), \(n\), 就得到 9 个方程: \[
\begin{pmatrix}
L_{1jkl}G_{1k,lj} + \delta(x) & L_{1jkl}G_{2k,lj} &
L_{1jkl}G_{3k,lj} \\
L_{2jkl}G_{1k,lj} & L_{2jkl}G_{2k,lj} + \delta(x) &
L_{2jkl}G_{3k,lj} \\
L_{3jkl}G_{1k,lj} & L_{3jkl}G_{2k,lj} & L_{3jkl}G_{3k,lj} +
\delta(x)
\end{pmatrix}
= \mathbf{0}_{3\times 3}
\] 方程的每一列代表了空间中三个方向上作用在 \(x=0\) 处的单位力, Green 函数 \(G_{nk}(x)\) 表示在 \(n\) 方向的单位集中力作用下, 弹性方程 \(k\) 方向的位移响应. 如果 \(L_{ijkl}\) 具有主对称性, 那么矩阵 \(K\) \(N\)
是对称矩阵, 所以 \[
G_{nk} = G_{kn}
\] 所以, \(G_{nk}(x)\)
同样也可以解释为: 在 \(k\)
方向的单位集中力作用下, 弹性方程 \(n\)
方向的位移响应. 这恰好与互易原理的结论相同.
各项同性材料的 Green
函数表示
对于各向同性材料, 方程左端的系数矩阵为 \[
\hat K \triangleq \mu \xi^{T} \xi ~\mathbf{I} + (\lambda + \mu)\xi
\xi^{T}
= \mu |\xi|^2 \mathbf{I}
+ (\lambda + \mu)
\begin{pmatrix}
\xi_1 \xi_1 & \xi_1 \xi_2 & \xi_1 \xi_3 \\
\xi_2 \xi_1 & \xi_2 \xi_2 & \xi_2 \xi_3 \\
\xi_3 \xi_1 & \xi_3 \xi_2 & \xi_3 \xi_3
\end{pmatrix}
\] 式中, \(\mathbf{I}\) 为 3
阶单位矩阵, \(\xi = (\xi_1, \xi_2,
\xi_3)^{\top}\), \(|\xi|^2 = \xi_1^2 +
\xi_2^2 + \xi_3^2\). 对于形式为 \(I +
\alpha vv^T\) 的矩阵, 其逆矩阵的形式依旧为 \(I + \beta vv^T\), 因此矩阵 \(\hat{K}\) 的逆为 \[
\hat{K}^{-1} = \frac{1}{\mu} \frac{1}{|\xi|^2} \mathbf{I}
- \frac{\lambda + \mu}{\mu (\lambda + 2\mu)} \frac{1}{|\xi|^4} \xi
\xi^{\top}, \quad
\det(\hat{K}) = (\lambda + 2\mu)\mu^2 |\xi|^6
\] 各项同性材料的 Green 函数为 \[
G(x) = \frac{1}{(2\pi)^3} \frac{1}{\mu} \int_{\mathbb{R}^3}
\frac{1}{|\xi|^2} \mathbf{I} ~\mathrm{d} \xi
- \frac{1}{(2\pi)^3} \frac{\lambda + \mu}{\mu (\lambda + 2\mu)}
\int_{\mathbb{R}^3} \frac{1}{|\xi|^4} \xi \xi^{\top} \mathrm{d} \xi
\] 将积分分成两部分进行计算 \[
I_0(x) = \int_{\mathbb{R}^3} \frac{1}{|\xi|^2} e^{i\xi \cdot x}
~\mathrm{d} \xi, \quad
I_{ij}(x) = - \frac{\partial^2}{\partial x_i \partial x_j}
\int_{\mathbb{R}^3} \frac{1}{|\xi|^4} e^{i\xi \cdot x} \mathrm{d} \xi
\] 一波操作得到 \[
I_0(x) = 2\pi^2 \frac{1}{|x|}, \quad
I_{ij} = \pi^2 \frac{\partial^2 |x|}{\partial x_i \partial x_j}
\] 代入方程中得到 \[
G_{ij}(x) = \frac{1}{4\pi \mu} \frac{1}{|x|} \delta_{ij}
- \frac{\lambda + \mu}{8\pi \mu (\lambda + 2\mu)}
\frac{\partial^2}{\partial x_i \partial x_j} |x|
\] > 对于径向对称幂函数 \(1/|\xi|^a\), Fourier
逆变换依旧得到径向对称幂函数 \(C_a/|x|^{d-a}\), 这时需要确定常数 \(C_a\), 当 \(a
< d\) 时, 有 > \[
> C_a=(2 \pi)^{\frac{n}{2}} \frac{2^{\frac{n-a}{2}}
\Gamma\left(\frac{n-a}{2}\right)}{2^{\frac{a}{2}}
\Gamma\left(\frac{a}{2}\right)} .
> \] > Radial
functions and the Fourier transform > > 对于 \(a > d\) 的情况, 积分是不收敛的, 但注意到
\(I_{ij}\) 还是关于 \(x_i\) 的二阶导数, 这就抵消了 2 幂次,
积分是收敛的. 尝试计算 \(I_{11}(x)\),
并令 \(x=(0,0,1)\), 最终得到 > \[
> C_{4} = 4\pi \int_{0}^{\infty} \frac{1}{r^3} \sin r - \frac{1}{r^2}
\cos r ~\mathrm{d}r
> \] > 有意思的是, 如果直接通过计算积分 \(\int_{\mathbb{R}^3} \frac{1}{|\xi|^4} e^{i\xi
\cdot x} \mathrm{d} \xi\) 求解常数 \(C_4\), 最终会得到上式积分项中关于 \(\sin r\) 部分 > \[
> C_4 = 4\pi \int_{0}^{\infty} \frac{1}{r^3} \sin r ~\mathrm{d}r
> \] > 而这个积分是不收敛的