Títol: Improved Numerical Methods for Elliptic Problems in Astrophysics

Autor/a: Adsuara Fuster, Jose Enrique

Resum: Les ePDEs (elliptic partial differential equations, en anglès) apareixen en una àmplia varietat d'àrees de les matemàtiques, la física i l'enginyeria. Són de particular interès en Astrofísica on apareixen, per exemple, quan es calcula el potencial gravitacional, en la solució de l'equació de Grad-Shafranov per magnetosferes lliures de forces, o d'imposar lligadures de divergència zero en la integració numèrica de les equacions MHD (magnetohydrodynamics, en anglès). En general, les ePDEs s'han de resoldre numèricament, establint una demanda cada vegada més gran d'algoritmes eficients i altament paral·lels per abordar la seua resolució computacional. El SRJ (scheduled relaxation Jacobi, en anglès) pertany a una prometedora classe de mètodes, atípic per la combinació de senzillesa i eficàcia, que s'ha introduït recentment per resoldre ePDEs lineals de tipus Poisson. És una extensió del mètode iteratiu clàssic de Jacobi utilitzat per resoldre sistemes d'equacions lineals del tipus Au = b. Hereta, d'entre altres, la seua robustesa. La seua metodologia es basa en el càlcul d'uns paràmetres apropiats per a una aproximació multinivell amb l'objectiu de minimitzar el nombre d'iteracions necessàries per a reduir el residual per davall d'una tolerància especificada. L'eficiència en la reducció del residual augmenta amb el nombre de nivells emprats en l'algoritme. Tanmateix, l'aplicació de la metodologia original per calcular els paràmetres d'estos esquemes SRJ òptims més enllà de 5 nivells és enormement dificultosa. Això és degut fonamentalment a la presència d'un sistema mixt algebraic-diferencial (no lineal) d'equacions el qual es torna cada vegada més rígid a mesura que augmenta el nombre de nivells. D'una banda, hem trobat una nova metodologia per a l'obtenció dels paràmetres dels esquemes òptims de l'algoritme SRJ que supera les limitacions de la metodologia original i proporciona els paràmetres per a estos esquemes amb un nombre elevat de nivells, fóra bo fins a 15, i per a resolucions de fins a 215 punts per dimensió. Això dóna lloc a factors d'acceleració de diversos centenars respecte del mètode de Jacobi en el cas de resolucions típiques i de milers en alguns casos amb altes resolucions. La major part de l'èxit en la recerca d'estos esquemes òptims amb més de 10 nivells es basa en una reducció analítica de la complexitat del sistema d'equacions abans esmentat. A més, s'estén l'algoritme original per aplicar-lo a certs sistemes d'equacions el·líptiques no lineals. D'altra banda, en un esquema típic SRJ, s'empra l'anterior conjunt de paràmetres en cicles de M iteracions consecutives fins que s'arriba a la tolerància prescrita. Presentem la forma analítica del conjunt òptim de factors de relaxació per al cas en què tots ells són estrictament diferents, i veiem que l'algoritme resultant és equivalent al mètode no estacionari de Richardson generalitzat, en el que es precondiciona la matriu del sistema d'equacions multiplicant per D = diag(A). El nostre mètode per estimar els pesos té l'avantatge que el càlcul explícit dels valors propis mínim i màxim de la matriu A (o la matriu d'iteració corresponent de l'esquema de Jacobi amb pes subjacent) es substitueix pel càlcul (molt més fàcil) de les freqüències mínima i màxima derivades de l'anàlisi d'estabilitat de von Neumann de l'operador el·líptic continu. Este conjunt de pesos també és l'òptim per al problema general, la qual cosa ens dóna la convergència més ràpida de tots els possibles esquemes SRJ per una estructura de malla donada. Ens referirem a ell com el mètode de Chebyshev-Jacobi. El factor d'amplificació del mètode es pot trobar analíticament i permet l'estimació exacta del nombre d'iteracions necessàries per a assolir la tolerància desitjada. També mostrem que a partir del conjunt de pesos calculats per l'esquema SRJ òptim per a una mida de cicle fix és possible calcular numèricament el valor òptim del factor de relaxació del mètode SOR (successive overrelaxation, en anglès) en alguns casos. Demostrem amb exemples pràctics, d'aplicació en Astrofísica, que el nostre mètode també funciona molt bé per als problemes de tipus Poisson en els que es fa servir una discretització d'alt ordre de l'operador Laplacià (per exemple, discretitzacions de 9- o 17- punts). Això té molt d'interès, ja que estes discretitzacions no produeixen matrius CO (consistently ordered, en anglès) i, per tant, la teoria de Young no es pot utilitzar per calcular el valor òptim del paràmetre de relaxació òptim de SOR. D'altra banda, els esquemes SRJ òptims deduïts ací són avantatjoses respecte a les implementacions existents per SOR pel que fa a discretitzacions d'alt ordre de l'operador Laplacià en la mesura que no cal recórrer als esquemes multicolors per a la seua execució en paral·lel. Presentem el mètode de Chebyshev-Jacobi fent servir una implementació purament MPI i una implementació híbrida OpenMP/MPI, ambdues sobre màquines de memòria compartida i de memòria distribuïda. Mostrem el seu rendiment i com escalen. També mostrem com arribar a velocitats de convergència notables amb execucions en paral·lel sobre GPUs quan la resolució d'equacions en derivades parcials el·líptiques amb diferències finites es fa utilitzant de manera conjunta el mètode de Chebyshev-Jacobi i les discretitzacions d'alt ordre. Finalment, tractar d'aplicar els nostres mètodes més enllà de l'àmbit de l'Astrofísica. En particular, abordem el problema de trobar els modes normals de vibració de l'ull humà. Este problema es pot resoldre amb una variant millorada de la metodologia que ací es presenta. La millora consisteix a estendre el càlcul del conjunt òptim de paràmetres al cas de matrius no definides positives. Les nostres idees sobre com procedir en este camp s'esbossen en el treball futur d'esta tesi.

URI: http://hdl.handle.net/10550/59545

Creado en: 2017

Relacionat Ítems a Google Scholar