notebook.community

Edit and run

Algoritmo de Crout

\begin{equation*} \begin{bmatrix} l_{11} & 0 & 0 & 0 \\ l_{21} & l_{22} & 0 & 0 \\ l_{31} & l_{32} & l_{33} & 0 \\ l_{41} & l_{42} & l_{43} & l_{44} \end{bmatrix} \begin{bmatrix} 1 & u_{12} & u_{13} & u_{14} \\ 0 & 1 & u_{23} & u_{24} \\ 0 & 0 & 1 & u_{34} \\ 0 & 0 & 0 & 1 \end{bmatrix} = \begin{bmatrix} a_{11} & a_{12} & a_{13} & a_{14} \\ a_{21} & a_{22} & a_{23} & a_{24} \\ a_{31} & a_{32} & a_{33} & a_{34} \\ a_{41} & a_{42} & a_{43} & a_{44} \end{bmatrix} \end{equation*}

Multiplicando

\begin{equation*} \begin{bmatrix} l_{11} & l_{11} u_{12} & l_{11} u_{13} & l_{11} u_{14} \\ l_{21} & l_{21} u_{12} + l_{22} & l_{21} u_{13} + l_{22} u_{23} & l_{21} u_{14} + l_{22} u_{24} \\ l_{31} & l_{31} u_{12} + l_{32} & l_{31} u_{13} + l_{32} u_{23} + l_{33} & l_{31} u_{14} + l_{32} u_{24} + l_{33} u_{34} \\ l_{41} & l_{41} u_{12} + l_{42} & l_{41} u_{13} + l_{42} u_{23} + l_{43} & l_{41} u_{14} + l_{42} u_{24} + l_{43} u_{34} + l_{44} \end{bmatrix} = \begin{bmatrix} a_{11} & a_{12} & a_{13} & a_{14} \\ a_{21} & a_{22} & a_{23} & a_{24} \\ a_{31} & a_{32} & a_{33} & a_{34} \\ a_{41} & a_{42} & a_{43} & a_{44} \end{bmatrix} \end{equation*}

Despejando y reemplazando

\begin{equation} \begin{bmatrix} a_{11} & 0 & 0 & 0 \\ a_{21} & a_{22} - l_{21} u_{12} & 0 & 0 \\ a_{31} & a_{32} - l_{31} u_{12} & a_{33} - l_{31} u_{13} - l_{32} u_{23} & 0 \\ a_{41} & a_{42} - l_{41} u_{12} & a_{43} - l_{41} u_{13} - l_{42} u_{23} & a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34} \end{bmatrix} \begin{bmatrix} 1 & \frac{a_{12}}{l_{11}} & \frac{a_{13}}{l_{11}} & \frac{a_{14}}{l_{11}} \\ 0 & 1 & \frac{a_{23} - l_{21} u_{13}}{l_{22}} & \frac{a_{24} - l_{21} u_{14}}{l_{22}} \\ 0 & 0 & 1 & \frac{a_{34} - l_{31} u_{14} - l_{32} u_{24}}{l_{33}} \\ 0 & 0 & 0 & 1 \end{bmatrix} = \begin{bmatrix} a_{11} & a_{12} & a_{13} & a_{14} \\ a_{21} & a_{22} & a_{23} & a_{24} \\ a_{31} & a_{32} & a_{33} & a_{34} \\ a_{41} & a_{42} & a_{43} & a_{44} \end{bmatrix} \end{equation}

Ejemplo

Factorizar la matriz $\mathbf{A}$

\begin{equation*} \begin{bmatrix} 1 & 1 & 2 & 3 \\ 2 & 1 & -1 & 1 \\ 3 & -1 & -1 & 2 \\ -1 & 2 & 3 & -1 \end{bmatrix} \end{equation*}

Calculando la primera columna de $\mathbf{L}$ y la primera fila $\mathbf{U}$

\begin{align*} l_{11} &= a_{11} = 1 \\ l_{21} &= a_{21} = 2 \\ l_{31} &= a_{31} = 3 \\ l_{41} &= a_{41} = -1 \\ u_{12} &= \frac{a_{12}}{l_{11}} = \frac{1}{1} = 1 \\ u_{13} &= \frac{a_{13}}{l_{11}} = \frac{2}{1} = 2 \\ u_{14} &= \frac{a_{14}}{l_{11}} = \frac{3}{1} = 3 \end{align*}

Reemplazando

\begin{equation*} \begin{bmatrix} 1 & 0 & 0 & 0 \\ 2 & & 0 & 0 \\ 3 & & & 0 \\ -1 & & & \end{bmatrix} \begin{bmatrix} 1 & 1 & 2 & 3 \\ 0 & 1 & & \\ 0 & 0 & 1 & \\ 0 & 0 & 0 & 1 \end{bmatrix} \end{equation*}

Calculando la segunda columna de $\mathbf{L}$ y la segunda fila $\mathbf{U}$

\begin{align*} l_{22} &= a_{22} - l_{21} u_{12} = 1 - 2 (1) = -1 \\ l_{32} &= a_{32} - l_{31} u_{12} = -1 - 3 (1) = -4 \\ l_{42} &= a_{42} - l_{41} u_{12} = 2 - (-1) (1) = 3 \\ u_{23} &= \frac{a_{23} - l_{21} u_{13}}{l_{22}} = \frac{-1 - 2 (2)}{-1} = 5 \\ u_{24} &= \frac{a_{24} - l_{21} u_{14}}{l_{22}} = \frac{1 - 2 (3)}{-1} = 5 \end{align*}

Reemplazando

\begin{equation*} \begin{bmatrix} 1 & 0 & 0 & 0 \\ 2 & -1 & 0 & 0 \\ 3 & -4 & & 0 \\ -1 & 3 & & \end{bmatrix} \begin{bmatrix} 1 & 1 & 2 & 3 \\ 0 & 1 & 5 & 5 \\ 0 & 0 & 1 & \\ 0 & 0 & 0 & 1 \end{bmatrix} \end{equation*}

Calculando la tercera columna de $\mathbf{L}$ y la tercera fila $\mathbf{U}$

\begin{align*} l_{33} &= a_{33} - l_{31} u_{13} - l_{32} u_{23} = -1 - 3 (2) - (-4) (5) = 13 \\ l_{43} &= a_{43} - l_{41} u_{13} - l_{42} u_{23} = 3 - (-1) (2) - 3 (5) = -10 \\ u_{34} &= \frac{a_{34} - l_{31} u_{14} - l_{32} u_{24}}{l_{33}} = \frac{2 - 3 (3) - (-4) (5)}{13} = 1 \end{align*}

Reemplazando

\begin{equation*} \begin{bmatrix} 1 & 0 & 0 & 0 \\ 2 & -1 & 0 & 0 \\ 3 & -4 & 13 & 0 \\ -1 & 3 & -10 & \end{bmatrix} \begin{bmatrix} 1 & 1 & 2 & 3 \\ 0 & 1 & 5 & 5 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 0 & 1 \end{bmatrix} \end{equation*}

Calculando la cuarta columna de $\mathbf{L}$

\begin{equation*} l_{44} = a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34} = -1 - (-1) (3) - 3 (5) - (-10) (1) = -3 \end{equation*}

Reemplazando

\begin{equation*} \begin{bmatrix} 1 & 0 & 0 & 0 \\ 2 & -1 & 0 & 0 \\ 3 & -4 & 13 & 0 \\ -1 & 3 & -10 & -3 \end{bmatrix} \begin{bmatrix} 1 & 1 & 2 & 3 \\ 0 & 1 & 5 & 5 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 0 & 1 \end{bmatrix} \end{equation*}

Patrón de cálculo

Primer patrón

\begin{equation*} \begin{array}{llll} \color{blue}{l_{11} = a_{11}} & u_{12} = \cfrac{a_{12}}{l_{11}} & u_{13} = \cfrac{a_{13}}{l_{11}} & u_{14} = \cfrac{a_{14}}{l_{11}} \\ \color{blue}{l_{21} = a_{21}} & l_{22} = a_{22} - l_{21} u_{12} & u_{23} = \cfrac{a_{23} - l_{21} u_{13}}{l_{22}} & u_{24} = \cfrac{a_{24} - l_{21} u_{14}}{l_{22}} \\ \color{blue}{l_{31} = a_{31}} & l_{32} = a_{32} - l_{31} u_{12} & l_{33} = a_{33} - l_{31} u_{13} - l_{32} u_{23} & u_{34} = \cfrac{a_{34} - l_{31} u_{14} - l_{32} u_{24}}{l_{33}} \\ \color{blue}{l_{41} = a_{41}} & l_{42} = a_{42} - l_{41} u_{12} & l_{43} = a_{43} - l_{41} u_{13} - l_{42} u_{23} & l_{44} = a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34} \end{array} \end{equation*}

Lo anterior puede escribirse como

\begin{equation*} l_{i1} = a_{i1} \end{equation*}

para $i=1, 2, 3, 4 = 1, \dots, m$

Segundo patrón

\begin{equation*} \begin{array}{llll} l_{11} = a_{11} & \color{blue}{u_{12} = \cfrac{a_{12}}{l_{11}}} & \color{blue}{u_{13} = \cfrac{a_{13}}{l_{11}}} & \color{blue}{u_{14} = \cfrac{a_{14}}{l_{11}}} \\ l_{21} = a_{21} & l_{22} = a_{22} - l_{21} u_{12} & u_{23} = \cfrac{a_{23} - l_{21} u_{13}}{l_{22}} & u_{24} = \cfrac{a_{24} - l_{21} u_{14}}{l_{22}} \\ l_{31} = a_{31} & l_{32} = a_{32} - l_{31} u_{12} & l_{33} = a_{33} - l_{31} u_{13} - l_{32} u_{23} & u_{34} = \cfrac{a_{34} - l_{31} u_{14} - l_{32} u_{24}}{l_{33}} \\ l_{41} = a_{41} & l_{42} = a_{42} - l_{41} u_{12} & l_{43} = a_{43} - l_{41} u_{13} - l_{42} u_{23} & l_{44} = a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34} \end{array} \end{equation*}

Lo anterior puede escribirse como

\begin{equation*} l_{1j} = \frac{a_{1j}}{l_{11}} \end{equation*}

para $j=2, 3, 4 = 2, \dots, n$

Tercer patrón

\begin{equation*} \begin{array}{llll} l_{11} = a_{11} & u_{12} = \cfrac{a_{12}}{l_{11}} & u_{13} = \cfrac{a_{13}}{l_{11}} & u_{14} = \cfrac{a_{14}}{l_{11}} \\ l_{21} = a_{21} & l_{\color{blue}{22}} = a_{\color{blue}{22}} - l_{\color{blue}{2}1} u_{1\color{blue}{2}} & u_{\color{blue}{2}3} = \cfrac{a_{\color{blue}{2}3} - l_{\color{blue}{2}1} u_{13}}{l_{\color{blue}{22}}} & u_{\color{blue}{2}4} = \cfrac{a_{\color{blue}{2}4} - l_{\color{blue}{2}1} u_{14}}{l_{\color{blue}{22}}} \\ l_{31} = a_{31} & l_{3\color{blue}{2}} = a_{3\color{blue}{2}} - l_{31} u_{1\color{blue}{2}} & l_{\color{green}{33}} = a_{\color{green}{33}} - l_{\color{green}{3}1} u_{1\color{green}{3}} - l_{\color{green}{3}2} u_{2\color{green}{3}} & u_{\color{green}{3}4} = \cfrac{a_{\color{green}{3}4} - l_{\color{green}{3}1} u_{14} - l_{\color{green}{3}2} u_{24}}{l_{\color{green}{33}}} \\ l_{41} = a_{41} & l_{4\color{blue}{2}} = a_{4\color{blue}{2}} - l_{41} u_{1\color{blue}{2}} & l_{4\color{green}{3}} = a_{4\color{green}{3}} - l_{41} u_{1\color{green}{3}} - l_{42} u_{2\color{green}{3}} & l_{44} = a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34} \end{array} \end{equation*}

Lo anterior puede escribirse como

\begin{align*} l_{jj} &= a_{jj} - \sum l_{j?} u_{?j} \\ l_{?j} &= a_{?j} - \sum l_{??} u_{?j} \\ u_{j?} &= \frac{a_{j?} - \sum l_{j?} u_{??}}{l_{jj}} \end{align*}

para $j=2, 3 = 2, \dots, n-1$

Cuarto patrón

\begin{equation*} \begin{array}{llll} l_{11} = a_{11} & u_{12} = \cfrac{a_{12}}{l_{11}} & u_{13} = \cfrac{a_{13}}{l_{11}} & u_{14} = \cfrac{a_{14}}{l_{11}} \\ l_{21} = a_{21} & l_{22} = a_{22} - l_{21} u_{12} & u_{2\color{blue}{3}} = \cfrac{a_{2\color{blue}{3}} - l_{21} u_{1\color{blue}{3}}}{l_{22}} & u_{2\color{blue}{4}} = \cfrac{a_{2\color{blue}{4}} - l_{21} u_{1\color{blue}{4}}}{l_{22}} \\ l_{31} = a_{31} & l_{\color{blue}{3}2} = a_{\color{blue}{3}2} - l_{\color{blue}{3}1} u_{12} & l_{33} = a_{33} - l_{31} u_{13} - l_{32} u_{23} & u_{3\color{green}{4}} = \cfrac{a_{3\color{green}{4}} - l_{31} u_{1\color{green}{4}} - l_{32} u_{2\color{green}{4}}}{l_{33}} \\ l_{41} = a_{41} & l_{\color{blue}{4}2} = a_{\color{blue}{4}2} - l_{\color{blue}{4}1} u_{12} & l_{\color{green}{4}3} = a_{\color{green}{4}3} - l_{\color{green}{4}1} u_{13} - l_{\color{green}{4}2} u_{23} & l_{44} = a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34} \end{array} \end{equation*}

Lo anterior puede escribirse como

\begin{align*} l_{ij} &= a_{ij} - \sum l_{i?} u_{?j} \\ u_{ji} &= \frac{a_{ji} - \sum l_{j?} u_{?i}}{l_{jj}} \end{align*}

para

\begin{align*} i &= 3, 4 \\ &= 4 \end{align*}

Quinto patrón

\begin{equation*} \begin{array}{llll} l_{11} = a_{11} & u_{12} = \cfrac{a_{12}}{l_{11}} & u_{13} = \cfrac{a_{13}}{l_{11}} & u_{14} = \cfrac{a_{14}}{l_{11}} \\ l_{21} = a_{21} & l_{22} = a_{22} - l_{2\color{blue}{1}} u_{\color{blue}{1}2} & u_{23} = \cfrac{a_{23} - l_{2\color{blue}{1}} u_{\color{blue}{1}3}}{l_{22}} & u_{24} = \cfrac{a_{24} - l_{2\color{blue}{1}} u_{\color{blue}{1}4}}{l_{22}} \\ l_{31} = a_{31} & l_{32} = a_{32} - l_{3\color{blue}{1}} u_{\color{blue}{1}2} & l_{33} = a_{33} - l_{3\color{green}{1}} u_{\color{green}{1}3} - l_{3\color{green}{2}} u_{\color{green}{2}3} & u_{34} = \cfrac{a_{34} - l_{3\color{green}{1}} u_{\color{green}{1}4} - l_{3\color{green}{2}} u_{\color{green}{2}4}}{l_{33}} \\ l_{41} = a_{41} & l_{42} = a_{42} - l_{4\color{blue}{1}} u_{\color{blue}{1}2} & l_{43} = a_{43} - l_{4\color{green}{1}} u_{\color{green}{1}3} - l_{4\color{green}{2}} u_{\color{green}{2}3} & l_{44} = a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34} \end{array} \end{equation*}

Lo anterior puede escribirse como

\begin{align*} l_{jj} &= a_{jj} - \sum l_{jk} u_{kj} \\ l_{ij} &= a_{ij} - \sum l_{ik} u_{kj} \\ u_{ji} &= \frac{a_{ji} - \sum l_{jk} u_{ki}}{l_{jj}} \end{align*}

para

\begin{align*} k &= 1 \\ &= 1, 2 \end{align*}

Sexto patrón

\begin{equation*} \begin{array}{llll} l_{11} = a_{11} & u_{12} = \cfrac{a_{12}}{l_{11}} & u_{13} = \cfrac{a_{13}}{l_{11}} & u_{14} = \cfrac{a_{14}}{l_{11}} \\ l_{21} = a_{21} & l_{22} = a_{22} - l_{21} u_{12} & u_{23} = \cfrac{a_{23} - l_{21} u_{13}}{l_{22}} & u_{24} = \cfrac{a_{24} - l_{21} u_{14}}{l_{22}} \\ l_{31} = a_{31} & l_{32} = a_{32} - l_{31} u_{12} & l_{33} = a_{33} - l_{31} u_{13} - l_{32} u_{23} & u_{34} = \cfrac{a_{34} - l_{31} u_{14} - l_{32} u_{24}}{l_{33}} \\ l_{41} = a_{41} & l_{42} = a_{42} - l_{41} u_{12} & l_{43} = a_{43} - l_{41} u_{13} - l_{42} u_{23} & \color{blue}{l_{44} = a_{44} - l_{41} u_{14} - l_{42} u_{24} - l_{43} u_{34}} \end{array} \end{equation*}

Lo anterior puede escribirse como

\begin{equation*} l_{nn} = a_{nn} - \sum_{k=1}^{n-1} l_{nk} u_{kn} \end{equation*}

Fórmula matemática

\begin{align*} i &= 1, \dots, m \\ & \quad l_{i1} = a_{i1} \\ j &= 2, \dots, n \\ & \quad u_{1j} = \frac{a_{1j}}{l_{11}} \\ j &= 2, \dots, n-1 \\ & \quad l_{jj} = a_{jj} - \sum_{k=1}^{j-1} l_{jk} u_{kj} \\ & \quad i = j+1, m \\ & \quad \quad l_{ij} = a_{ij} - \sum_{k=1}^{j-1} l_{ik} u_{kj} \\ & \quad \quad u_{ji} = \frac{a_{ji} - \sum_{k=1}^{j-1} l_{jk} u_{ki}}{l_{jj}} \\ l_{nn} &= a_{nn} - \sum_{k=1}^{n-1} l_{nk} u_{kn} \end{align*}

Seudocódigo

function lu_crout(A)
    for i=1 to m do
        l(i,1) = a(i,1)
    end for
    for j=2 to n do
        u(1,j) = a(1,j)/l(1,1)
    end for
    for j=2 to n-1 do
        suma = a(j,j)
        for k=1 to j-1 do
            suma = suma - l(j,k)*u(k,j)
        end for
        l(j,j) = suma
        for i=j+1 to m do
            sumav = a(i,j)
            sumah = a(j,i)
            for k=1 to j-1 do
                sumav = sumav - l(i,k)*u(k,j)
                sumah = sumah - l(j,k)*u(k,i)
            end for
            l(i,j) = sumav
            u(j,i) = sumah/l(j,j)
        end for
    end for
    suma = a(n,n)
    for k=1 to n-1 do
        suma = suma - l(n,k)*u(k,n)
    end for
    l(n,n) = suma
    return a
end function

Implementación descomposición de Crout





In [5]:



    


import numpy as np

def lu_crout(a):
    m, n = a.shape
    u = np.eye(m)
    l = np.zeros((m,n))
    for i in range(0,m):
        l[i,0] = a[i,0]
                 
    for j in range(1,n):
        u[0,j] = a[0,j]/l[0,0]
        
    for j in range(1,n-1):
        suma = a[j,j]
        for k in range(0,j):
            suma = suma - l[j,k]*u[k,j]
        l[j,j] = suma
        for i in range(j+1,m):
            sumav = a[i,j]
            sumah = a[j,i]
            for k in range(0,j):
                sumav = sumav - l[i,k]*u[k,j]
                sumah = sumah - l[j,k]*u[k,i]
            l[i,j] = sumav
            u[j,i] = sumah/l[j,j]
            
    suma = a[n-1,n-1]
    for k in range(0,n-1):
        suma = suma - l[n-1,k]*u[k,n-1]
    l[n-1,n-1] = suma
    
    return l, u



    











In [6]:



    


A = np.array([[1,1,2,3],
              [2,1,-1,1],
              [3,-1,-1,2],
              [-1,2,3,-1]],float)
print(A)



    


















    






[[ 1.  1.  2.  3.]
 [ 2.  1. -1.  1.]
 [ 3. -1. -1.  2.]
 [-1.  2.  3. -1.]]

















In [7]:



    


#resultado
L, U = lu_crout(A)
print(L)
print(U)



    


















    






[[  1.   0.   0.   0.]
 [  2.  -1.   0.   0.]
 [  3.  -4.  13.   0.]
 [ -1.   3. -10.  -3.]]
[[ 1.  1.  2.  3.]
 [ 0.  1.  5.  5.]
 [ 0.  0.  1.  1.]
 [ 0.  0.  0.  1.]]

















In [8]:



    


#revisando el resultado
np.dot(L,U)



    


















    
Out[8]:








array([[ 1.,  1.,  2.,  3.],
       [ 2.,  1., -1.,  1.],
       [ 3., -1., -1.,  2.],
       [-1.,  2.,  3., -1.]])

















In [9]:



    


A = np.array([[3,-0.1,-0.2],
              [0.1,7,-0.3],
              [0.3,-0.2,10]],float)
print(A)



    


















    






[[  3.   -0.1  -0.2]
 [  0.1   7.   -0.3]
 [  0.3  -0.2  10. ]]

















In [10]:



    


#resultado
L, U = lu_crout(A)
print(L)
print(U)



    


















    






[[  3.           0.           0.        ]
 [  0.1          7.00333333   0.        ]
 [  0.3         -0.19        10.01204188]]
[[ 1.         -0.03333333 -0.06666667]
 [ 0.          1.         -0.04188482]
 [ 0.          0.          1.        ]]

















In [11]:



    


#revisando el resultado
np.dot(L,U)



    


















    
Out[11]:








array([[  3. ,  -0.1,  -0.2],
       [  0.1,   7. ,  -0.3],
       [  0.3,  -0.2,  10. ]])

Content source: ClaudioVZ/Metodos_numericos_I

Similar notebooks:

notebook.community | gallery | about