def lapack_plu ( ):

## lapack_plu() demonstrates tril() and triu() by uncompressing LAPACK output.
#
  import numpy as np

  print ( '' )
  print ( 'lapack_plu():' )
  print ( '  Demonstrate the use of tril() and triu() functions.' )
  print ( '  Reconstruct P, L, U factors from compressed LAPACK output.' )
  print ( '' )
  print ( '  Given the following matrix A:' )
  A = np.array ( [ [ 1, 2, 3 ], [ 4, 5, 6 ], [ 7, 8, 0 ] ] )
  n = A.shape[0]
  print ( A )
  print ( "  Calling LAPACK for a trans(P) * L * U factorization" )
  print ( '    call dgetrf ( a, n, n, lda, ipiv, info )' )
  print ( '  returns the following compressed information:' )
  print ( '' )
  print ( '    compressed LU info:' )
  ALU = np.array ( [ \
    [ 7.0,      8.0,      0.0 ], \
    [ 0.142857, 0.857143, 3.0 ], \
    [ 0.571429, 0.500000, 4.50000 ] ] )
  print ( ALU )
  print ( '' )
  print ( '    compressed pivot info P:' )
  ipiv = np.array ( [ 3, 3, 3 ], dtype = int )
  print ( ipiv )
#
#  Reconstruct P, L, U matrices.
#
  print ( '' )
  print ( '  U = np.triu ( ALU, 0 )' )
  U = np.triu ( ALU, 0 )
  print ( U )
  print ( '' )
  print ( '  L = np.tril ( ALU, -1 ) + np.identity ( n )' )
  L = np.tril ( ALU, -1 ) + np.identity ( n )
  print ( L )
  print ( '' )
  print ( '  P is a permutation, reconstructed from p' )
  P = np.identity ( n )
  for i in range ( 0, n ):
    k = ipiv[i] - 1
    if ( k != i ):
      P[[k,i]] = P[[i,k]]
  print ( P )
#
#  Verify factorization.
#
  print ( '' )
  print ( "  Check PLU = trans(P) * L * U = A?" )
  PLU = np.matmul ( np.transpose ( P ), np.matmul ( L, U ) )
  print ( PLU )

  return

if ( __name__ == "__main__" ):
  lapack_plu ( )