Mathematical Engineering - Matematica Numerica

Exercise 06

Divided by topic

MATEMATICA NUMERICA A.A. 2018 - 2019 Ingegneria Matematica Prof. A. Quarteroni Prof. A. Manzoni, Dr. I. Fumagalli Esercitazione 6 Metodi diretti per sistemi lineari (2) Esercizio 1 Si consideri la matrice A = " b 1c conb; c2R, e0< "1, ovvero molto piccolo. a)Si applichi il metodo di eliminazione di Gauss (senza pivoting) per determinare la fattorizzazione LUdiA. b)Si applichi il metodo di eliminazione di Gauss con pivotazione per righe per determinare la fattoriz-zazione LU diA. c)Consideriamo la matrice non singolare A =2 41 1 + 0 :5
10 15 3 2 2 20 3 6 43 5; si esegua con Matlab la fattorizzazione LU mediante il metodo di eliminazione di Gauss senza pivoting per righe (usando la funzionelugauss.m) e con pivoting per righe (usando il comando Matlablu). Quale dierenza si nota ? Sapreste motivarla ? Esercizio 2 Utilizzare opportunamente i comandi Matlabonesediagper costruire una matriceA2Rn n conn= 20 tale che,a 1;i= 1 ea i;1= 1 peri= 1; : : : ; n,a ii= 4 peri= 2; : : : ; ned innea i;i+1= 1ea i+1;i= 1per i= 2; : : : ; n1. Tutti gli altri coecienti diAsono nulli. Tale matrice si dicesparsapoiché il numero degli elementi non-nulli diAèO(n)n2 , dunque molto elevato1 . Per visualizzare la posizione degli elementi non nulli diAsi utilizzi il comandospy. Calcolare inoltre la fattorizzazioneLUdiA. Utilizzare il comandospyper determinare se le matriciLed Usiano o meno sparse. A vostro parere, le caratteristiche di sparsità dei fattoriLeUsono vantaggiose o svantaggiose ?1. In tal caso è vantaggioso memorizzare solo gli elementi non-nulli di A. In Matlab ciò viene fatto attraverso il comando sparse, ad esempioB=sparse(A).1 Esercizio 3 Per limitare il fenomeno del riempimento (l l-in) è possibile eettuare operazioni preliminari di permutazione tra righe, che prendono il nome direordering2 . L'algoritmo di Cuthill-McKee (e la sua variante inversa) permette di permutare una matrice sparsa simmetrica in una matrice a banda più compatta. La variante inversa di tale algoritmo è implementata in Matlab nel comandosymrcm; in particolare, r = symrcm(S) restituisce una permutazionertale cheS(r,r)tende ad avere tutti gli elementi non nulli vicini alla diagonale. a)Generare con il comandoB = buckyuna matriceB2R60 60 dallagallerydi Matlab. Vericare che tale matrice ha una banda piuttosto larga. b)Eettuare la fattorizzazione LU della matriceB, e mostrare il pattern dei fattori L ed U. c)Eettuare il reordering della matriceBcon il comandosymrcme, indicando conBrla matrice riordinata, determinare la fattorizzazione LU di tale matrice. Mostrare il pattern dei fattori L ed U, confrontando i risultati con quelli ottenuti al punto precedente. Esercizio 4 Si consideri una matriceAtridiagonale (caso particolare di matrice a banda) e invertibile : A =2 6 6 6 6 4a 1c 10 e2a 2: :: :::c n1 0e na n3 7 7 7 7 5; dovefa ign i=1, fc ign 1 i=1, fe ign i=22 R. Si vuole fattorizzare tale matrice nel prodottoA = LU(caso particolare difattorizzazione LU) doveLeUsono due matrici bidiagonali, della forma L =2 6 6 6 41 0 21 :::: :: 0 n13 7 7 7 5; U =2 6 6 6 6 4 1 10 2: :: ::: n1 0 n3 7 7 7 7 5; allo scopo di risolvere il sistema lineareAx=bin modo eciente. Essendo infattiLeUcasi particolari di matrici rispettivamente triangolari inferiori e triangolari superiori, tale riscrittura permette di semplicare la risoluzione del sistema, adottando gli algoritmi disostituzione in avantiesostituzione al l'indietro. Questi algoritmi prevedono la risoluzione sequenziale di ciascuna riga di un sistema triangolarenn, andando rispettivamente, dalla prima all'ultima, o dall'ultima alla prima. a)Fornire l'espressione dei coecienti i,  ie iin funzione dei coecienti di A. b)Applicare le formule trovate in a) e costruire un algoritmo per risolvere il sistema lineareAx=b, con b= [b 1; : : : ; b n]T 2Rn , mediante sostituzioni in avanti e all'indietro (algoritmo di Thomas). c)Quante operazionioating-pointsono richieste dall'algoritmo di Thomas ? d)Implementare in Matlab, mediante unafunction, l'algoritmo di Thomas per risolvere un generico sistema lineare tridiagonale. L'intestazione dellafunctiondovrà essere : function [L,U,x] = thomas(A,b)2. Da non confondersi con il pivoting2 Esercizio 5 Sia data la matrice A =2 47 1 3 1 8 2 3 2 93 5; ed il vettoreb= Ax exdove x ex= [1 ;1;1]T . a)Mostrare che la matriceAammette un'unica fattorizzazioneLUsenza ricorrere alla pivotazione. b)Calcolare, tramite il comandoludi Matlab, la fattorizzazioneLUdiA. Questo comando esegue la pivotazione e quindi restituisce le matriciL,UePtali cheLU = PA. Vericare che in questo caso specico Matlab non ha eseguito alcuna permutazione. Giusticare la risposta alla luce della teoria. c)Utilizzare gli algoritmi di sostituzione in avanti ed all'indietro, implementati nelle functionsfwsube bksubrispettivamente, per risolvere il sistema lineareAx=battraverso la fattorizzazioneLU. d)Mostrare che la matriceAammette un'unica fattorizzazione di Cholesky,A = HT HdoveHè una matrice triangolare superiore ; trovareHtramite il comando Matlabchol. Utilizzare gli algoritmi di sostituzione in avanti ed all'indietro, per risolvere il sistema lineareAx=battraverso la fattorizzazione di Cholesky. Esercizio 6 Si consideri la seguente matrice A =2 41 0 2 1 0 1 33 5: a)Fornire un intervallo di valori ditale che la matriceAammetta fattorizzazione di Cholesky. b)Calcolare analiticamente la fattorizzazione di Cholesky di A in funzione dinell'insieme dei valori ammissibili trovato al punto precendente. Vericare il calcolo nel caso= 1utilizzando il comando choldi Matlab. c)Sempre nel caso= 1, risolvere in Matlab i tre sistemi lineari seguenti, utilizzando opportunamente la fattorizzazione di Cholesky trovata al punto precedente e le funzionifwsubebksub: Ax i= e i; i = 1;2;3; dovee 1= [1 ;0;0]T ,e 2= [0 ;1;0]T ee 3= [0 ;0;1]T . Utilizzare questi risultati per calcolareA 1 . Si spieghi perché la fattorizzazione di Cholesky è conveniente dal punto di vista computazionale per il calcolo dell'inversa diA.3 Esercizio 7 Si consideri un sistema sico formato daN+1molle ideali, disposte in la su un piano orizzontale e agganciate l'una all'altra tramite degli anelli ; la prima e l'ultima molla sono ancorate ad una parete. La distanza fra le due pareti èL, e tutte le molle hanno la stessa costante di rigidezzaK. Si può mostrare che la congurazione di equilibrio delle molle è ottenuta risolvendo il sistema lineare K2 6 6 6 6 6 4 2 1 12 1 :::: ::: :: 12 1 123 7 7 7 7 7 52 6 6 6 6 6 6 4x 1 : : : : : : xN1 xN3 7 7 7 7 7 7 5=2 6 6 6 6 6 6 40 : : : : : : 0 K L3 7 7 7 7 7 7 5 dovex irappresenta la posizione dell' i-esimo anello, coni= 1; : : : ; N;x 0= X N+1= 0 poiché la prima e l'ultima molla sono ancorate ad una parete. Confrontare i tempi di calcolo dei metodi\,LU, Thomas, Cholesky sul sistema lineare associato alla matric sistema delle molle, facendo variare il valore del parametronche denisce la dimensione della matriceAda n= 200an= 2000, con passo200. Per quale valore dinle dierenze sono apprezzabili ?4