using Plots

# nice viz for matrices
function lookat(A; redrow=0, rounding=2, showtext=true)
  n = size(A,1)
  plot(legend=false, axis=false)
  rowcolor = redrow > 0 ? :red : :black  
  for i=1:n, j=1:n  
      scatter!( [j],[i], ann=  showtext ? (j,i,round(A[i,j],digits=rounding), :white ) : (j,i,"") ,
                color=abs(A[i,j]) > .0001 ? (i==redrow ? rowcolor : :black) : :white, 
                marker=:square, markersize=30, aspectratio=1, yflip=true, yaxis=[.5,n+.5],xaxis=[.5,n+.5])
  end
  plot!()
end

A = [2 12 4 7; 29 12 69 24 ; 7 7 7 7; 3 7 5 2]
A = A*1.0



lookat(A)

L = fill(0.0, 4,4)

L[2,1] = A[2,1]/A[1,1]

A[2,:] -= L[2,1]*A[1,:]; # subtract the multiplier times the first row from the second row
lookat(A, redrow=2)

L[3,1] = A[3,1]/A[1,1] # the diagonal entry that we divide by is called the "pivot"

A[3,:] -= L[3,1]*A[1,:]; # subtract the multiplier times the first row from the second row
lookat(A, redrow=3)

L[4,1] = A[4,1]/A[1,1] 
A[4,:] -= L[4,1]*A[1,:]; # subtract the multiplier times the first row from the second row
lookat(A, redrow=4)

L[3,2] = A[3,2]/A[2,2] 
A[3,:] -= L[3,2]*A[2,:]; # subtract the multiplier times the first row from the second row
lookat(A, redrow=3)

L[4,2] = A[4,2]/A[2,2] 
A[4,:] -= L[4,2]*A[2,:]; # subtract the multiplier times the first row from the second row
lookat(A, redrow=4)

L[4,3] = A[4,3]/A[3,3] 
A[4,:] -= L[4,3]*A[3,:]; # subtract the multiplier times the first row from the second row
lookat(A, redrow=4)

n = 6
A = rand(n,n)
keepA = [copy(A)]
row = [0]
L = fill(0.0,n,n)
for j=1:n, i=j+1:n
    L[i,j] = A[i,j] / A[j,j]
    A[i,:] -= L[i,j] * A[j,:]
    push!(keepA,copy(A))
    push!(row,i)
end

using Interact

@manipulate for i =1: length(keepA)
    lookat( keepA[i], redrow=row[i])
end

using LinearAlgebra

A = keepA[1]

U = keepA[end]

L

L += I

A

L * U

A ≈ L * U

A = [2 12 4 7; 29 12 69 24 ; 7 7 7 7; 3 7 5 2]

A'


A= [ 1 2 3; 4 5 6]

A'



dump(A')

A = rand( 1:9, 2,3,4)

permutedims(A, [3 2 1])

inv( rand(4,6))

A = rand(3,3)

B = rand(3,3)

inv( A*B )

inv(B) * inv(A)

x = 1:10


x * x'

x =[ 1 ,2 ,4, 5]

y = [ 2, 4, 6]

x * y'