Matrix factorizations and direct solution of linear systers 51-9 12.GV96,p.115j Algorithm 3:GEPP/PLU decomposition of a rectangular matrix(outer product) InDt:AG设m Output:L E Rmxm(unit lower triangular matrix) PeRgmpw (P is represented with an index vector p such that y=Pz=) L←lm:U←0∈Rmx";p=1,2,3,,m:andr←1； fork=1 to n Findμsuch that≤≤n and A=maxrsism Ak ir.the change An.kin Ar,tn,Lpir-1 Lr.tr-1,and pe pr r+l:m.r t Ar+im.k/Ark for i=n |A←A-LU r4r+1 Algorithm 4:GEPP/PLU decomposition of a rectangular matrix(gaxpy) nput:A∈Rm Output:L ERxm(unit lower triangular matrix). UERx"(upper triangular matrix-row echelon form),and PER (permutation matrix)so that PA=LU (Pis represented with an index vector that records row interchanges In∈row r and row之rwere interchanged in step】 L- ；andr←-1； for j=I to n hange vv元 Solve the triangular system,L1.1=vir- Update v ←Vrm-Lrm.1r-1z Find u such that=maxrsismlvl ifu≠0，then Exchange var for i=l tor-1. Exchange←+Lrd Lr+1:m.r←Vr+lm/r ←r+1 Matrix Factorizations and Direct Solution of Linear Systems 51-9 12. [GV96, p. 115] Algorithm 3: GEPP/PLU decomposition of a rectangular matrix (outer product) Input: A ∈ R m×n Output: L ∈ R m×m (unit lower triangular matrix) U ∈ R m×n (upper triangular matrix - row echelon form) P ∈ R m×m (permutation matrix) so that P A = LU (P is represented with an index vector p such that y = Pz ⇔ yj = zpj ) L ← Im; U ← 0 ∈ R m×n ; p = [1, 2, 3, . . . , m]; and r ← 1; for k = 1 to n Find µ such that r ≤ µ ≤ m and |Aµ,k| = maxr≤i≤m |Ai,k| if Aµ,k 6= 0, then Exchange Aµ,k:n ↔ Ar,k:n, Lµ,1:r−1 ↔ Lr,1:r−1, and pµ ↔ pr Lr+1:m,r ← Ar+1:m,k/Ar,k Ur,k:n ← Ar,k:n for i = r + 1 to m for j = k + 1 to n Ai,j ← Ai,j − Li,rUr,j r ← r + 1 Algorithm 4: GEPP/PLU decomposition of a rectangular matrix (gaxpy) Input: A ∈ R m×n Output: L ∈ R m×m (unit lower triangular matrix), U ∈ R m×n (upper triangular matrix - row echelon form), and P ∈ R m×m (permutation matrix) so that P A = LU (P is represented with an index vector π that records row interchanges πr = µ means row r and row µ ≥ r were interchanged in step r) L ← Im ∈ R m×m; U ← 0 ∈ R m×n ; and r ← 1; for j = 1 to n v ← A1:m,j if r > 1, then for i = 1 to r − 1, Exchange vi ↔ vπi Solve the triangular system, L1:r−1,1:r−1 · z = v1:r−1 U1:r−1,j ← z Update vr:m ← vr:m − Lr:m,1:r−1 · z Find µ such that |vµ| = maxr≤i≤m |vi| if vµ 6= 0, then πr ← µ Exchange vµ ↔ vr for i = 1 to r − 1, Exchange Lµ,i ↔ Lr,i Lr+1:m,r ← vr+1:m/vr Ur,j ← vr r ← r + 1