Mètode multipol ràpid

El mètode multipol ràpid (FMM) és una tècnica numèrica que es va desenvolupar per accelerar el càlcul de forces de llarg abast en el problema de n cossos. Això ho fa ampliant la funció de Green del sistema mitjançant una expansió multipolar, que permet agrupar fonts que es troben a prop i tractar-les com si fossin una única font.

El FMM també s'ha aplicat per accelerar el solucionador iteratiu en el mètode dels moments (MOM) tal com s'aplica als problemes d'electromagnètica computacional, i en particular en el bioelectromagnetisme computacional. El FMM va ser introduït per primera vegada d'aquesta manera per Leslie Greengard i Vladimir Rokhlin Jr. i es basa en l'expansió multipolar de l'equació vectorial de Helmholtz. En tractar les interaccions entre funcions de base llunyanes mitjançant el FMM, els elements de la matriu corresponents no cal que s'emmagatzemen explícitament, el que resulta en una reducció significativa de la memòria requerida. Si el FMM s'aplica de manera jeràrquica, pot millorar la complexitat dels productes de vectors matrius en un solucionador iteratiu de ${\mathcal {O}}(N^{2})$ a ${\mathcal {O}}(N)$ en aritmètica finita, és a dir, donada una tolerància $\varepsilon$, es garanteix que el producte matriu-vector es troba dins d'una tolerància $\varepsilon .$ La dependència de la complexitat de la tolerància $\varepsilon$ és ${\mathcal {O}}(\log(1/\varepsilon ))$, és a dir, la complexitat de FMM és ${\mathcal {O}}(N\log(1/\varepsilon ))$. Això ha ampliat l'àrea d'aplicabilitat del MOM a problemes molt més grans dels que eren possibles anteriorment.

S'ha dit que el FMM, introduït per Rokhlin Jr. i Greengard, és un dels deu millors algorismes del segle XX. L'algorisme FMM redueix la complexitat de la multiplicació matriu-vector que implica un determinat tipus de matriu densa que pot sorgir de molts sistemes físics.

El FMM també s'ha aplicat per tractar de manera eficient la interacció de Coulomb en el mètode Hartree-Fock i els càlculs de la teoria funcional de la densitat en química quàntica.

Mètode multipol ràpid: interpolació d'un pol a x=3 amb un polinomi de Chebyshev d'ordre 5

Esbós de l'algoritme

En la seva forma més simple, el mètode multipol ràpid busca avaluar la funció següent:

$f(y)=\sum _{\alpha =1}^{N}{\frac {\phi _{\alpha }}{y-x_{\alpha }}}$

on $x_{\alpha }\in [-1,1]$ són un conjunt de pols i $\phi _{\alpha }\in \mathbb {C}$ són els pesos dels pols corresponents en un conjunt de punts $\{y_{1},\ldots ,y_{M}\}$ amb $y_{\beta }\in [-1,1]$. Aquesta és la forma unidimensional del problema, però l'algoritme es pot generalitzar fàcilment a múltiples dimensions i nuclis diferents de $(y-x)^{-1}$.

Ingènuament, avaluant $f(y)$ activat $M$ punts requereix ${\mathcal {O}}(MN)$ operacions. L'observació crucial darrere del mètode multipol ràpid és que si la distància entre $y$ i $x$ és prou gran, doncs $(y-x)^{-1}$ està ben aproximat per un polinomi. Concretament, deixa $-1<t_{1}<\ldots <t_{p}<1$ ser els nodes de l'ordre de Txebixev $p\geq 2$ i deixar $u_{1}(y),\ldots ,u_{p}(y)$ siguin els polinomis de base de Lagrange corresponents. Es pot demostrar que el polinomi interpolador:

${\frac {1}{y-x}}=\sum _{i=1}^{p}{\frac {1}{t_{i}-x}}u_{i}(y)+\epsilon _{p}(y)$

convergeix ràpidament amb l'ordre polinomial, $|\epsilon _{p(y)}|<5^{-p}$, sempre que el pol estigui prou lluny de la regió d'interpolació, $|x|\geq 3$ i $|y|<1$. Això es coneix com a "expansió local".

L'acceleració del mètode multipolar ràpid deriva d'aquesta interpolació: sempre que tots els pols estiguin "allunyats", avaluem la suma només als nodes de Chebyshev a un cost de ${\mathcal {O}}(Np)$, i després interpolar-lo a tots els punts desitjats a un cost de ${\mathcal {O}}(Mp)$ :

$\sum _{\alpha =1}^{N}{\frac {\phi _{\alpha }}{y_{\beta }-x_{\alpha }}}=\sum _{i=1}^{p}u_{i}(y_{\beta })\sum _{\alpha =1}^{N}{\frac {1}{t_{i}-x_{\alpha }}}\phi _{\alpha }$

Des que $p=-\log _{5}\epsilon$, on $\epsilon$ és la tolerància numèrica, el cost total és ${\mathcal {O}}((M+N)\log(1/\epsilon ))$.

Per garantir que els pols estiguin ben separats, es subdivideix recursivament l'interval unitari de manera que només ${\mathcal {O}}(p)$ els pols acaben en cada interval. A continuació, s'utilitza la fórmula explícita dins de cada interval i la interpolació per a tots els intervals que estan ben separats. Això no fa malbé l'escala, ja que un necessita com a màxim $\log(1/\epsilon )$ nivells dins de la tolerància donada.