﻿ 边界元算法在计算地球动力学中的应用
 中国科学院大学学报  2016, Vol. 33 Issue (1): 89-96 PDF

Applications of boundary-element method in computational geodynamics
LI Zhonghai
State Key Laboratory of Continental Tectonics and Dynamics, Institute of Geology, Chinese Academy of Geological Sciences, Beijing 100037, China
Abstract: The boundary-element method (BEM) is widely used in the engineering and technology fields. However,applications of BEM in geosciences are relatively rare at the present. In this paper, we first describe the boundary-integral equations in Stokes flow problems with the Green's functions under several boundary conditions. Then, we discuss applications of BEM numerical modeling in plate subduction dynamics and subduction-zone seismic anisotropy. The results indicate that BEM can be used for solving many difficult problems in geodynamics.
Key words: boundary-element method     boundary-integral equation     Green' s function     geodynamics

1 斯托克斯方程的格林函数 1.1 无限空间的格林函数

 ${{\partial }_{i}}{{\sigma }_{ij}}={{F}_{j}}\delta (x{{x}_{0}}),\partial {{}_{i}}{{u}_{i}}=0$ (1)

 ${{u}_{i}}={{J}_{ij}}{{F}_{j}}/\eta ,{{\sigma }_{ik}}={{K}_{ijk}}{{F}_{}}$ (2)

 ${{\partial }_{i}}{{J}_{ij}}=0,{{}_{i}}{{K}_{ijk}}={{\delta }_{jk}}\delta (x{{x}_{0}}).$ (3)

 ${{J}_{ij}}\left( r \right)=\frac{1}{8\pi }\left( \frac{{{\delta }_{ij}}}{r}+\frac{{{r}_{i}}{{r}_{j}}}{{{r}^{3}}} \right),$ (4)
 ${{P}_{j}}\left( r \right)=\frac{3}{4\pi }\frac{{{r}_{j}}}{{{r}^{3}}},$ (5)
 ${{K}_{ijk}}\left( r \right)=\frac{1}{4\pi }\frac{{{r}_{i}}{{r}_{j}}{{r}_{k}}}{{{r}^{5}}}$ (6)

1.2 一个自由滑动边界的格林函数

 \begin{align} & {{G}_{ij}}(x-{{x}_{0}})={{{\hat{G}}}_{ij}}(x-{{x}_{0}})+G_{ij}^{IM}(x-x_{0}^{IM})= \\ & {{{\hat{G}}}_{ij}}(x-{{x}_{0}})+R_{pj}^{\left( 3 \right)}G_{ip}^{IM}(x-x_{0}^{IM}) \\ \end{align} (7)

 ${{R}^{\left( 3 \right)}}_{ij}=\left( \begin{matrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & -1 \\ \end{matrix} \right).$ (8)
1.3 一个刚性边界的格林函数

 \begin{align} & {{J}_{ij}}={{{\hat{J}}}_{ij}}+{{{\overset{\lower0.5em\hbox{\smash{\scriptscriptstyle\smile}$}}{J}}}_{ij}}, \\ & {{P}_{j}}={{{\hat{P}}}_{^{j}}}+{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{P}}}_{j}}, \\ & {{K}_{ijk}}={{{\hat{K}}}_{ijk}}+{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}}}{K}}}_{ijk}}, \\ \end{align} (9)

 \begin{align} & {{{\overset{\lower0.5em\hbox{\smash{\scriptscriptstyle\frown}}}{J}}}_{ij}}=\frac{1}{8\pi }\left\{ \frac{{{\delta }_{ij}}}{R}-\frac{{{R}_{i}}{{R}_{j}}}{{{R}^{3}}}+2{{x}_{03}}{{\Delta }_{j}}\left[ {{x}_{03}}\left( \frac{{{\delta }_{ij}}}{{{R}^{3}}}-\frac{3{{R}_{i}}{{R}_{j}}}{{{R}^{5}}} \right)+\frac{{{\delta }_{i3}}}{{{R}^{3}}}{{R}_{j}} \right] \right. \\ & -\frac{1}{{{R}^{3}}}({{\delta }_{ij}}{{R}_{3}}+{{\delta }_{3j}}{{R}_{i}})+\left. \frac{3{{R}_{i}}{{R}_{j}}{{R}_{k}}}{{{R}^{5}}} \right\}, \\ \end{align} (10)
 ${{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{P}}}_{j}}=\frac{1}{4\pi }-\left[ \frac{{{R}_{j}}}{{{R}^{3}}}-2{{x}_{03}}{{\Delta }_{j}}\left( \frac{{{\delta }_{3j}}}{{{R}^{3}}}-\frac{3{{R}_{j}}{{R}_{3}}}{{{R}^{5}}} \right) \right],$ (11)
 \begin{align} & {{{\overset{\lower0.5em\hbox{\smash{\scriptscriptstyle\frown}}}{K}}}_{ijk}}=-\frac{3}{4\pi }\left\{ -\frac{{{R}_{i}}{{R}_{j}}{{R}_{k}}}{{{R}^{5}}} \right.2{{x}_{03}}{{\Delta }_{j}}\left[ -\frac{{{x}_{03}}}{{{R}^{5}}} \right.{{\delta }_{ik}}{{R}_{j}}+\frac{{{x}_{3}}}{{{R}^{5}}}({{\delta }_{jk}}{{R}_{i}}+{{\delta }_{ij}}{{R}_{k}}) \\ & +\frac{{{R}_{i}}{{\delta }_{3j}}{{R}_{k}}}{{{R}^{5}}}\left. \left. \frac{{{R}_{i}}{{R}_{j}}{{R}_{k}}{{x}_{3}}}{{{R}^{7}}} \right] \right\} \\ \end{align} (12)

1.4 两个平行的刚性边界的格林函数

 ${{J}_{ij}}={{v}_{ij}}+{{w}_{ij}}$ (13)
 ${{P}_{j}}={{q}_{j}}+{{s}_{j}}$ (14)

1.4.1 速度和压强格林函数的基本项

 \begin{align} & {{v}_{ij}}=\frac{1}{4\pi }\left[ {{\delta }_{ij}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \lambda \rho \right){{B}_{1}}\lambda d\lambda +{{\delta }_{i\alpha }}{{\delta }_{j\beta }}{{r}_{\alpha }}{{r}_{\beta }}\frac{1}{\rho } \right.{{\int }^{\infty }}_{0} \\ & \lambda {{J}_{1}}\lambda \rho {{B}_{2}}\lambda d\lambda {{\delta }_{i3}}{{\delta }_{j3}}{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\lambda \rho {{B}_{3}}\lambda d\lambda +sgn \\ & \left( {{x}_{3}}h \right){{\delta }_{i3}}{{\delta }_{j\alpha }}+{{\delta }_{i\alpha }}{{\delta }_{j3}}{{r}_{\alpha }}{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\lambda \rho {{B}_{4}}\lambda d\left. \lambda \right] \\ \end{align} (15)
 \begin{align} & {{q}_{j}}=\frac{1}{2\pi }\left[ {{\delta }_{j\alpha }}{{r}_{\alpha }}\frac{1}{\rho } \right.{{\int }^{\infty }}_{0}\lambda {{J}_{1}}\lambda \rho {{B}_{1}}\lambda d\lambda +sgn\left( {{x}_{3}}h \right){{\delta }_{j3}} \\ & {{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{4}}\lambda d\left. \left( \lambda \right) \right] \\ \end{align} (16)

1) 当x3≥h时,满足

 \begin{align} & {{B}_{1}}\left( \lambda \right)=\frac{sinh\lambda h}{sinh\lambda H}sinh\lambda (H{{x}_{3}}), \\ & {{B}_{2}}\left( \lambda \right)={{B}_{1}}\left( \lambda \right), \\ & {{B}_{3}}\left( \lambda \right)=\frac{d}{d\lambda }\left[ \frac{sinh\lambda h}{sinh\lambda H} \right.sinh\lambda (H{{x}_{3}}), \\ & {{B}_{4}}\left( \lambda \right)=\left. \frac{sinh\lambda h}{sinh\lambda H}cosh\lambda (H{{x}_{3}}) \right]. \\ \end{align} (17)

2) 当x3<h时,h需要替换为(H－h), x3替换为(H－x3), 可得

 \begin{align} & {{B}_{1}}\left( \lambda \right)=\frac{sinh\lambda \left( Hh \right)}{sinh\lambda H}sinh\lambda {{x}_{3}}, \\ & {{B}_{2}}\left( \lambda \right)={{B}_{1}}\left( \lambda \right), \\ & {{B}_{3}}\left( \lambda \right)=\frac{d}{d\lambda }\left[ \frac{sinh\lambda \left( Hh \right)}{sinh\lambda H}sinh\lambda {{x}_{3}} \right], \\ & {{B}_{4}}\left( \lambda \right)=\frac{sinh\lambda \left( Hh \right)}{sinh\lambda H}cosh\lambda {{x}_{3}}. \\ \end{align} (18)
1.4.2 速度和压强格林函数的补充项

 \begin{align} & {{w}_{\alpha \beta }}=\frac{1}{4\pi }\frac{\partial }{\partial {{r}_{\beta }}}\frac{{{r}_{\alpha }}}{\rho }{{\int }^{\infty }}_{0}\xi {{J}_{1}}\left( \rho \xi \right){{A}_{1}}\left( \xi \right)d\xi = \\ & -\frac{{{\delta }_{\alpha \beta }}}{4\pi }\frac{1}{\rho }{{\int }^{\infty }}_{0}\xi {{J}_{1}}\left( \rho \xi \right){{A}_{1}}\left( \xi \right)d\xi -\frac{{{r}_{\alpha }}{{r}_{\beta }}}{2\pi }\frac{1}{\partial {{\rho }^{2}}}\frac{1}{\rho }{{\int }^{\infty }}_{0}\xi {{J}_{1}}\left( \rho \xi \right){{A}_{1}}\left( \xi \right)d\xi \\ \end{align} (19)

 \begin{align} & {{A}_{1}}\left( \xi \right)=\frac{1}{sin{{h}^{2}}\left( \xi H \right){{\left( \xi H \right)}^{2}}}\{\xi hH{{x}_{3}}sinh\xi \left( H{{x}_{3}} \right)sinh\xi \left( Hh \right)+ \\ & {{x}_{3}}[hsinh\xi Hcosh\xi (H-{{x}_{3}}-h)-Hsinh\xi hcosh\xi {{x}_{3}}]+ \\ & \xi H{{x}_{3}}sinh\xi Hcosh\xi (H{{x}_{3}})\frac{d}{d\xi }\left( \frac{sinh\xi \left( Hh \right)}{sinh\xi H} \right) \\ & Hsinh\xi {{x}_{3}}\left[ sinh\xi H\frac{d}{d\xi } \right]\left( \frac{sinh\xi h}{sinh\xi H} \right)+\xi H\frac{d}{d\xi }\left( \frac{sinh\xi \left( Hh \right)}{sinh\xi H} \right). \\ \end{align} (20)
 ${{w}_{33}}=\frac{1}{4\pi }{{\int }^{\infty }}_{0}{{J}_{0}}\left( \rho \xi \right){{A}_{2}}\xi d\xi .$ (21)
 \begin{align} & {{A}_{2}}\left( \xi \right)=\frac{{{\xi }^{2}}}{sin{{h}^{2}}\left( \xi H \right){{\left( \xi H \right)}^{2}}}\{{{x}_{3}}[Hsinh\xi hcosh\xi {{x}_{3}}hsinh\xi Hcosh\xi \\ & (H{{x}_{3}}h)+\xi hHsinh\xi (H{{x}_{3}})sinh\xi \left( Hh \right)+\xi Hcosh\xi (H{{x}_{3}}) \\ & sinh\xi H\frac{d}{d\xi }\left( \frac{sinh\xi \left( Hh \right)}{sinh\xi H} \right)\xi {{H}^{2}}sinh\xi {{x}_{3}}\frac{d}{d\xi }\left( \frac{sinh\xi \left( Hh \right)}{sinh\xi H} \right)+ \\ & Hsinh\xi {{x}_{3}}sinh\xi H\frac{d}{d\xi }\left( \frac{sinh\xi h}{sinh\xi H} \right). \\ \end{align} (22)
 \begin{align} & {{w}_{3\alpha }}{{_{(}}_{\alpha 3}}_{)}=-\frac{{{r}_{\alpha }}}{2\pi }\left[ \frac{\partial }{\partial {{\rho }^{2}}} \right.{{\int }^{\infty }}_{0}{{J}_{0}}\left( \rho \xi \right){{A}_{3}}\left( \xi \right)d\xi + \\ & \left. ({{\delta }_{i3}}{{\delta }_{j3}})\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \rho \xi \right){{A}_{4}}\left( \xi \right)d\xi \right] \\ \end{align} (23)
 \begin{align} & {{A}_{3}}\xi =\frac{1}{sin{{h}^{2}}\left( \xi H \right)}\left[ hsinh\xi Hcosh\xi \right.(H{{x}_{3}}h)Hsinh\left( \xi h \right)cosh\left( \xi {{x}_{3}} \right)+ \\ & \frac{\xi }{sin{{h}^{2}}\xi h\xi {{H}^{2}}}\left\{ {{x}_{3}}\xi H \right.hsinh\xi {{x}_{3}}h+\frac{Hsinh\xi hsinh\xi \left( H-{{x}_{3}} \right)}{sinh\xi H} \\ & -\text{ }\frac{\xi h{{H}^{2}}sinh\xi {{x}_{3}}sinh\xi \left( Hh \right)}{sinh\xi H} \\ \end{align} (24)

 \begin{align} & {{A}_{4}}\left( \xi \right)=\frac{\xi }{sin{{h}^{2}}\left( \xi H \right){{\left( \xi H \right)}^{2}}}-H\left( H-h \right)sinh\xi hsinh\xi {{x}_{3}}+ \\ & {{x}_{3}}hsinh\xi Hsinh\xi (H{{x}_{3}}h)+{{x}_{3}}Hsinh\xi hsinh\xi {{x}_{3}}] \\ \end{align} (25)
 ${{s}_{j}}=\frac{1}{2\pi }\left[ {{\delta }_{j\alpha }}{{r}_{\alpha }}\frac{1}{\rho }{{\int }^{\infty }}_{0}{{J}_{1}}\left( \rho \xi \right){{A}_{5}}\left( \xi \right)d\xi +{{\delta }_{j3}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \rho \xi \right){{A}_{6}}\left( \xi \right)d\xi \right].$ (26)

 \begin{align} & {{A}_{5}}\left( \xi \right)=\frac{{{\xi }^{2}}}{sin{{h}^{2}}\left( \xi H \right){{\left( \xi H \right)}^{2}}} \\ & \left\{ sinh\xi Hsinh\xi (H{{x}_{3}})\left[ sinh\xi H\frac{d}{d\xi } \right. \right.\left( \frac{sinh\xi h}{sinh\xi H} \right)+ \\ & \xi H\frac{d}{d\xi }\left( \frac{sinh\xi \left( Hh \right)}{sinh\xi H} \right)+cosh\xi (H{{x}_{3}}) \\ & [\xi hHsinh\xi \left( Hh \right)+\left( Hh \right)sinh\xi hsinh\xi H]\}, \\ \end{align} (27)
 \begin{align} & {{A}_{6}}\left( \xi \right)=\frac{{{\xi }^{2}}}{sin{{h}^{2}}\left( \xi H \right){{\left( \xi H \right)}^{2}}} \\ & \left\{ sinh\xi Hcosh\xi (H{{x}_{3}})\left[ sinh\xi H \right. \right.\frac{d}{d\xi }\left( \frac{sinh\xi h}{sinh\xi H} \right) \\ & \xi H\frac{d}{d\xi }\left( \frac{sinh\xi \left( Hh \right)}{sinh\xi H} \right)+sinh\xi (H{{x}_{3}}) \\ & [\xi hHsinh\xi \left( Hh \right)\left( Hh \right)sinh\xi hsinh\xi H]\}. \\ \end{align} (28)
1.4.3 应力的格林函数

 ${{K}_{ijk}}={{P}_{j}}{{\delta }_{ik}}+\frac{\partial {{J}_{ij}}}{\partial {{x}_{k}}}+\frac{\partial {{J}_{kj}}}{\partial {{x}_{i}}}.$ (29)

 \begin{align} & {{K}_{\alpha \beta \gamma }}=\frac{{{\delta }_{\alpha \gamma }}{{r}_{\beta }}}{2\pi }\left[ \frac{1}{\rho } \right.{{\int }^{\infty }}_{0}\lambda {{J}_{1}}\left( \lambda \rho \right){{B}_{2}}\left( \lambda \right)d\lambda \frac{1}{\rho }{{\int }^{\infty }}_{0}\lambda {{J}_{1}}\left( \lambda \rho \right){{B}_{1}}\left( \lambda \right)d\lambda \\ & \frac{1}{\rho }{{\int }^{\infty }}_{0}{{J}_{1}}\left( \xi \rho \right){{A}_{5}}\left( \xi \right)d\xi 2\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}\xi {{J}_{1}}\left( \xi \rho \right){{A}_{1}}\left( \xi \right)\left. d\xi \right]\frac{{{\delta }_{\beta \gamma }}{{r}_{\alpha }}+{{\delta }_{\alpha \beta }}{{r}_{\gamma }}}{4\pi } \\ & \left[ \frac{1}{\rho } \right.{{\int }^{\infty }}_{0}\lambda {{J}_{1}}\left( \lambda \rho \right){{B}_{2}}\left( \lambda \right)d\lambda +2\frac{{{\partial }^{2}}}{{{\left( \partial {{\rho }^{2}} \right)}^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \lambda \rho \right){{B}_{1}}\left( \lambda \right)d\lambda \\ & 4\frac{\partial }{\partial {{\rho }^{2}}}\frac{1}{\rho }{{\int }^{\infty }}_{0}\xi {{J}_{1}}\xi \rho {{A}_{1}}\xi d\xi +{{r}_{\alpha }}{{r}_{\beta }}{{r}_{\gamma }}\pi \\ & \frac{\partial }{\partial {{\rho }^{2}}}\frac{1}{\rho }{{\int }^{\infty }}_{0}\lambda {{J}_{1}}\lambda \rho {{B}_{2}}\lambda d\lambda 2\frac{{{\partial }^{2}}}{{{(\partial {{\rho }^{2}})}^{2}}~}\left. \frac{1}{\rho }{{\int }^{\infty }}_{0}\xi {{J}_{1}}\left( \xi \rho \right){{A}_{1}}\left( \xi \right)d\xi \right] \\ \end{align} (30)
 \begin{align} & {{K}_{\alpha \beta 3}}=\frac{{{\delta }_{\alpha \beta }}}{4\pi }\left[ \frac{\partial }{\partial {{x}_{3}}} \right.{{\int }^{\infty }}_{0}{{J}_{0}}\left( \lambda \rho \right){{B}_{1}}\left( \lambda \right)d\lambda \\ & \frac{\partial }{\partial {{x}_{3}}}\frac{1}{\rho }{{\int }^{\infty }}_{0}\xi {{J}_{1}}\left( \xi \rho \right){{A}_{1}}\left( \xi \right)d\xi +{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{4}}\left( \lambda \right)d\lambda \\ & -2\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right)({{A}_{3}}\left( \xi \right)+{{A}_{4}}\left( \xi \right))d\xi \\ & +\frac{{{r}_{\alpha }}{{r}_{\beta }}}{4\pi }\frac{\partial }{\partial {{x}_{3}}}\frac{1}{\rho }{{\int }^{\infty }}_{0}\lambda {{J}_{1}}\left( \lambda \rho \right){{B}_{2}}\left( \lambda \right)d\lambda - \\ & 2\frac{\partial }{\partial {{x}_{3}}\partial {{\rho }^{2}}}\frac{1}{\rho }{{\int }^{\infty }}_{0}\xi {{J}_{1}}\left( \xi \rho \right){{A}_{1}}\left( \xi \right)d\xi + \\ & 2\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{4}}\left( \lambda \right)d\lambda 4\frac{{{\partial }^{2}}}{{{\left( \partial {{\rho }^{2}} \right)}^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right) \\ & \left. ({{A}_{3}}\left( \xi \right)+{{A}_{4}}\left( \xi \right))d\xi \right] \\ \end{align} (31)
 \begin{align} & {{K}_{3\beta 3}}=\frac{{{r}_{\beta }}}{2\pi }\left[ \frac{\partial }{\partial {{x}_{3}}} \right]{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{4}}\left( \lambda \right)d\lambda \\ & 2\frac{{{\partial }^{2}}}{\partial {{x}_{3}}{{\rho }^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right)({{A}_{3}}\left( \xi \right)+{{A}_{4}}\left( \xi \right))d\xi \\ & \left. \frac{1}{\rho }{{\int }^{\infty }}_{0}\lambda {{J}_{1}}\left( \lambda \rho \right){{B}_{1}}\left( \lambda \right)d\lambda -1\rho {{\int }^{\infty }}_{0}{{J}_{1}}\left( \xi \rho \right){{A}_{5}}\left( \xi \right)d\xi \right] \\ \end{align} (32)
 \begin{align} & {{K}_{\alpha 3\gamma }}=\frac{{{r}_{\alpha }}{{r}_{\gamma }}}{\pi }\left[ \frac{\partial }{\partial {{\rho }^{2}}} \right]{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{4}}\left( \lambda \right)d\lambda - \\ & \left. 2\frac{{{\partial }^{2}}}{{{\left( \partial {{\rho }^{2}} \right)}^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right)({{A}_{3}}\left( \xi \right){{A}_{4}}\left( \xi \right))d\xi \right]+\frac{{{\delta }_{\alpha \gamma }}}{2\pi } \\ & [2\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right)({{A}_{3}}\left( \xi \right){{A}_{4}}\left( \xi \right))d\xi {{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right){{A}_{6}}\left( \xi \right)d\xi ], \\ \end{align} (33)
 \begin{align} & {{K}_{\alpha 33}}=\frac{{{r}_{\alpha }}}{4\pi }\frac{\partial }{\partial {{x}_{3}}}{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{4}}\left( \lambda \right)d\lambda + \\ & 2\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right){{A}_{2}}\left( \xi \right)d\xi 2\frac{{{\partial }^{2}}}{\partial {{x}_{3}}\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right)({{A}_{3}}\left( \xi \right){{A}_{4}}\left( \xi \right))d\xi + \\ & 2\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \lambda \rho \right){{B}_{1}}\left( \lambda \right)d\lambda 2\frac{\partial }{\partial {{\rho }^{2}}}{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{3}}\left( \lambda \right)d\lambda , \\ \end{align} (34)
 \begin{align} & {{K}_{333}}=\frac{1}{2\pi }-{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{4}}\left( \lambda \right)d\lambda - \\ & {{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right){{A}_{6}}\left( \xi \right)d\xi +\frac{\partial }{\partial {{x}_{3}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \xi \rho \right){{A}_{2}}\left( \xi \right)d\xi + \\ & \left. \frac{\partial }{\partial {{x}_{3}}}{{\int }^{\infty }}_{0}{{J}_{0}}\left( \lambda \rho \right){{B}_{1}}\left( \lambda \right)d\lambda -\frac{\partial }{\partial {{x}_{3}}}{{\int }^{\infty }}_{0}\lambda {{J}_{0}}\left( \lambda \rho \right){{B}_{3}}\left( \lambda \right)d\lambda \right]. \\ \end{align} (35)
1.5 一个自由滑动和一个刚性边界的格林函数

 ${{G}_{ij}}(x{{x}_{0}})={{G}^{2w}}_{ij}(x{{x}_{0}})+{{G}^{IM}}_{ij}(x{{x}^{IM}}_{0})$ (36)

1.6 格林函数复杂项的预计算及数据库构建

2 边界积分方程

 \begin{align} & {{\chi }_{1}}({{x}_{0}}){{u}^{(1)}}_{j}({{x}_{0}})+\gamma {{\chi }_{2}}({{x}_{0}}){{u}^{(2)}}_{j}({{x}_{0}})= \\ & \frac{{{g}_{k}}\delta \rho }{{{\eta }_{1}}}\int x{{\prime }_{k}}{{J}_{ij}}(x{{x}_{0}}){{n}_{i}}\left( x \right)dS\left( x \right)+\left( 1\gamma \right) \\ & \int {{u}_{i}}\left( x \right){{K}_{ijk}}(x{{x}_{0}}){{n}_{k}}\left( x \right)dS\left( x \right) \\ \end{align} (37)

2.1 目标点在曲面S

 \begin{align} & \frac{1}{2}\left( 1+\gamma \right){{u}_{j}}({{x}_{0}})= \\ & -\frac{{{g}_{k}}\delta \rho }{{{\eta }_{1}}}\int x{{\prime }_{k}}{{J}_{ij}}(x{{x}_{0}}){{n}_{i}}\left( x \right) \\ & dS\left( x \right)+\left( 1\gamma \right)\int {{u}_{i}}\left( x \right){{K}_{ijk}}(x{{x}_{0}}){{n}_{k}}\left( x \right)dS\left( x \right). \\ \end{align} (38)

 \begin{align} & {{u}_{j}}({{x}_{0}})={{P}_{j}}({{x}_{0}})+{{Q}_{j}}({{x}_{0}}),{{P}_{j}}({{x}_{0}})= \\ & \int (z{{z}_{0}}){{J}_{ij}}(x{{x}_{0}}){{n}_{i}}\left( x \right)dS\left( x \right), \\ & {{Q}_{j}}({{x}_{0}})=\left( 1\gamma \right)\int [{{u}_{i}}\left( x \right){{u}_{i}}({{x}_{0}})]{{K}_{ijk}}(x{{x}_{0}}){{n}_{k}}\left( x \right)dS\left( x \right) \\ \end{align} (39)

2.2 目标点在流体V1或V2

 \begin{align} & {{u}_{j}}({{x}_{0}})={{P}_{j}}({{x}_{0}})+{{Q}_{j}}({{x}_{0}}),{{P}_{j}}({{x}_{0}})= \\ & 1\gamma \int (z{{z}_{0}}){{J}_{ij}}(x{{x}_{0}}){{n}_{i}}\left( x \right)dS\left( x \right),{{Q}_{j}}({{x}_{0}})= \\ & 1\gamma \gamma \int {{u}_{i}}\left( x \right){{K}_{ijk}}(x{{x}_{0}}){{n}_{k}}\left( x \right)dS\left( x \right). \\ \end{align} (40)

3 方程求解及误差分析

 Download: JPG larger image 图 1 (a)球体表面的曲面三角形网格化;(b)误差与网格单元数的关系 Fig. 1 (a) Surface discretization of a sphere using curved triangles; (b) dependence of the RMS error on the number of cells
4 应用实例: 板块俯冲带动力学

 Download: JPG larger image 图 2 板块俯冲模式依赖于板块与周围地幔的黏滞系数比值 Fig. 2 Dependence of subduction modes on the viscosity ratio between plate and surrounding mantle

 Download: JPG larger image 图 3 俯冲带地幔流动、物质变形及地震波各向异性的数值模拟原理示意图 Fig. 3 Numerical modeling methods for subduction-induced mantel flow,material deformation,and seismic anisotropy

5 小结

 [1] Jaswon M A. Integral equation methods in potential theory, Part 1[J]. Proceedings of the Royal Society of London Series A , 1963, 275 :23–32. DOI:10.1098/rspa.1963.0152 [2] Symm G T. Integral equation methods in potential theory, Part 2[J]. Proceedings of Royal Society of London Series A , 1963, 275 :33–46. DOI:10.1098/rspa.1963.0153 [3] Rizzo F J. An integral equation approach to boundary value problems of classical elastostatics[J]. Quarterly of Applied Mathematics , 1967, 25 :83–95. [4] Cruse T A, Rizzo F J. A direct formulation and numerical solution of the general transient elastodynamic problem, part 1[J]. Journal of Mathematical Analysis and Applications , 1968, 22 (1) :244–259. DOI:10.1016/0022-247X(68)90171-6 [5] Cruse T A, Rizzo F J. A direct formulation and numerical solution of the general transient elastodynamic problem, Part 2[J]. Journal of Mathematical Analysis and Applications , 1968, 22 (2) :341–355. DOI:10.1016/0022-247X(68)90177-7 [6] Brebbia C A. The boundary element method for engineers[M]. London: Pentech Press, 1978 : 189 . [7] Pozrikidis C. A practical guide to boundary element methods with the software library BEMLIB[M]. Boca Raton: CRC Press, 2002 : 440 . [8] Rjasanow S, Steinbach O. The fast solution of boundary integral equations[M]. New York: Springer US, 2007 : 284 . [9] Liu Y. Fast multipole boundary element method-Theory and applications in engineering[M]. Cambridge: Cambridge University Press, 2009 : 254 . [10] Li Z H, Ribe N M. Dynamics of free subduction from 3-D boundary-element modeling[J]. Journal of Geophysical Research: Solid Earth , 2012, 117 . [11] Li Z H, Di Leo J F, Ribe N M. Subduction-induced mantle flow, finite strain and seismic anisotropy: Numerical modeling[J]. Journal of Geophysical Research: Solid Earth , 2014, 119 :5052–5076. DOI:10.1002/2014JB010996 [12] Pozrikidis C. Boundary integral and singularity methods for linearized viscous flow[M]. New York: Cambridge University Press, 1992 : 259 . [13] Blake J R. A note on the image system for a stokeslet in a no-slip boundary[J]. Mathematical Proceedings of the Cambridge Philosophical Society , 1971, 70 :303–310. DOI:10.1017/S0305004100049902 [14] Ascoli E P, Dandy D S, Leal L G. Buoyancy-driven motion of a deformable drop toward a planar wall at low Reynolds number[J]. Journal of Fluid Mechanics , 1990, 213 :287–311. DOI:10.1017/S0022112090002336 [15] Liron N, Mochon S. Stokes flow for a stokeslet between two parallel flat plates[J]. Journal of Engineering Mathematics , 1976, 10 :287–303. DOI:10.1007/BF01535565 [16] Pozrikidis C. The deformation of a liquid drop moving normal to a plane wall[J]. Journal of Fluid Mechanics , 1990, 215 :331–363. DOI:10.1017/S0022112090002671 [17] Manga M, Stone H A. Buoyancy-driven interaction between two deformable viscous drops[J]. Journal of Fluid Mechanics , 1993, 256 :647–683. DOI:10.1017/S0022112093002915 [18] Fraysse V, Giraud L, Gratton S, et al. A set of GMRES routines for real and complex arithmetics on high performance computers . CERFACS Technical Report TR/PA/03/3, 2001: 20. http://cn.bing.com/academic/profile?id=2080710049&encoded=0&v=paper_preview&mkt=zh-cn [19] Batchelor G K. An introduction to fluid dynamics[M]. New York: Cambridge University Press, 1967 : 615 . [20] Brenner H. The slow motion of a sphere through a viscous fluid towards a plane surface[J]. Chemical Engineering Science , 1961, 16 :242–251. DOI:10.1016/0009-2509(61)80035-3