Algebra lineare numerica per il calcolo ad alte prestazioni

The heart of numerical applications from a very broad range of domains, from statistics to engineering, geophysics, medical imaging and many others, consists in the fast and accurate solution of large size linear systems in many flavors. This is a classic topic in numerical linear algebra of renewed interest if, like in the present work, we address to large ill-posed problems, in which the matrix is possibly rank-deficient and the solution may not be unique, or to emerging hardware technologies like many-cores architectures such as Graphic Processing Units (GPUs). The first part of this thesis is devoted to the efficient solution of numerically rank-deficient least squares and nonnegative least squares problems. Such problems are particularly challenging and their exact solution would require exhaustive-search algorithms and thus their complexity is combinatorial, with a consequent difficulty to achieve high performances in computations. To this aim, we devise a new block column selection strategy, we investigate the derived methods with respect to accuracy of the solution found from the theoretical point of view and analyze their performance against the state-of-the-art algorithms, showing results with an extensive campaign of experiments. In particular, the nonnegative least squares solver presented naturally produces sparse solutions and it can be employed in many applications in the field of compressed sensing. As a relevant application, we present a numerical package for the computation of compressed multivariate near G-optimal finite probability measures for polynomial regression design of general degree by means of nonnegative least squares problems. In the second part of this thesis, we target the challenging question of how linear solvers for sparse matrices, intended here as matrices mostly filled up with zeros, can be adapted to modern and future computing facilities by proposing numerical methods specifically designed for the highly parallel and hybrid hardware platforms available nowadays. We first address the efficient design and implementation of parallel direct solvers, that is based upon matrix factorization, for a particular class of linear systems with block structured matrices, arising e.g. in optimal control problems. In particular, we discuss the use of novel parallel programming paradigms on existing direct solvers in order to improve their scalability, intended as the ability to solve problems of larger and larger size. When dealing with sparse linear systems, iterative methods are generally appreciated to be more suited for parallel architectures, as they are based upon simple computations such as matrix-vector multiplication. These methods are often coupled with a preconditioner in order to improve the rate of convergence. One of the most popular and powerful general purpose is the ILU preconditioner, that relies on the direct solution of sparse triangular systems, an inherently sequential task, and its application turns out to be the most time consuming operation on parallel architectures. In this work we deal with the solution of the algebraic problems arising from the simulation of incompressible flows with variable density and we show how an adequate iterative approach to the parallel solution of the triangular systems that result from the ILU preconditioner, turns out to be robust and efficient on massively parallel architectures.

La prima parte di questa tesi è dedicata alla risoluzione efficiente di problemi ai minimi quadrati dove la matrice della funzione obiettivo può non avere rango pieno, e di problemi ai minimi quadrati con vincoli di nonnegatività. Queste classi di problemi sono particolarmente impegnative poiché la loro risoluzione esatta richiede tecniche di ricerca enumerativa, data la loro natura combitatorica, con conseguenti difficoltà nel derivare algoritmi ad alte prestazioni. A tal fine deriviamo un algoritmo a blocchi per la selezione di colonne, investighiamo le proprietà dei metodi ottenuti da un punto di vista teorico in termini di accuratezza della soluzione e analizziamo le performance ottenute rispetto agli algoritmi stato dell'arte, mostrando i risultati di un'esaustiva campagna di esperimenti. In particolare, il solutore ai minimi quadrati nonnegativi qui proposto favorisce naturalmente la sparsità della soluzione e perciò può essere adoperato in applicazioni di compressed sensing. Come applicazione rilavante, presentiamo un pacchetto numerico per la compressione quasi G-ottimale di misure di probabilità finite in più dimensioni per la regressione polinomiale di grado elevato basato sulla risoluzioni di problemi ai minimi quadrati nonnegativi. La seconda parte della tesi affronta l'impegnativa questione su come adattare i solutori per sistemi lineari con matrici sparse, intese qui come quelle matrici nelle quali la maggior parte degli elementi sono nulli, alle architetture moderne e future proponendo metodi specificatamente pensati per le piattaforme hardware ibride e massivamente parallele oggigiorno disponibili. Ci dedichiamo al design e all'implementazione efficiente di solutori lineari di tipo diretto, cioè basati su fattorizzazioni matriciali, per la risoluzione di una classe particolare di sistemi lineari con matrici a blocchi derivanti, ad esempio, dalla risoluzione numerica di problemi di controllo ottimo. Nello specifico, discutiamo l'utilizzo di paradigmi di programmazione parallela innovativi applicati ai solutori esistenti per migliorarne la scalabilità, cioè la capacità di gestire una mole di lavoro crescente. Per la risoluzione di sistemi lineari sparsi, i metodi iterativi sono generalmente ritenuti più adatti per l'implementazione su architetture parallele in quanto basati su operazioni semplici, quali il prodotto matrice-vettore. Questi metodi vengono spesso accoppiati con una tecnica di precondizionamento, per migliorarne la velocità di convergenza. Uno dei precondizionatori generici più diffusi ed efficaci è il precondizionatore ILU, basato sulla risoluzione di sistemi triagolari sparsi, un calcolo intrinsecamente sequentaziale, e la sua applicazione risulta essere l'operazione più dispendiosa sulle architetture parallele. In questo lavoro affrontiamo la risoluzione dei problemi algebrici derivanti dalla simulazione di flussi incomprimibili con densità variabile e mostriamo che un adeguato approccio iterativo per la risoluzione dei sistemi triangolari risultanti dall'applicazione del precondizionatore ILU si rivela robusto ed efficiente sulle architetture parallele.