Advanced Search

Question: About the robustness of Student:-NumericalAnalysis:-MatrixDecomposition

Question:About the robustness of Student:-NumericalAnalysis:-MatrixDecomposition

sursumCorda 1244 Maple 2025
linear-algebra factorization decomposition
July 31 2025
1

As the following worksheet shows, Student:-NumericalAnalysis:-MatrixDecomposition cannot factorize the input matrix m and throws an error, but if we simply reorder or exchange the elements of m, no error will be raised. (The reason for setting 'method' = 'LDLt' is that LinearAlgebra:-LUDecomposition can be used for other methods.) 
 

> 

restart
> 

with(Student:-NumericalAnalysis, MatrixDecomposition)
> 

m := Matrix([[3*(sqrt(3)+1)/8,-1/2,1/2,-(sqrt(3)+1)/8,-1/2,1/2,-(sqrt(3)+1)/8,-1/2,1/2,-(sqrt(3)+1)/8],
> 

             [-1/2,sqrt(3)-1,-(sqrt(3)-1),-1/2,0,0,1/2,0,0,1/2],
> 

             [1/2,-(sqrt(3)-1),sqrt(3)-1,1/2,0,0,-1/2,0,0,-1/2],
> 

             [-(sqrt(3)+1)/8,-1/2,1/2,3*(sqrt(3)+1)/8,1/2,-1/2,-(sqrt(3)+1)/8,1/2,-1/2,-(sqrt(3)+1)/8],
> 

             [-1/2,0,0,1/2,sqrt(3)-1,-(sqrt(3)-1),-1/2,0,0,1/2],
> 

             [1/2,0,0,-1/2,-(sqrt(3)-1),sqrt(3)-1,1/2,0,0,-1/2],
> 

             [-(sqrt(3)+1)/8,1/2,-1/2,-(sqrt(3)+1)/8,-1/2,1/2,3*(sqrt(3)+1)/8,1/2,-1/2,-(sqrt(3)+1)/8],
> 

             [-1/2,0,0,1/2,0,0,1/2,sqrt(3)-1,-(sqrt(3)-1),-1/2],
> 

             [1/2,0,0,-1/2,0,0,-1/2,-(sqrt(3)-1),sqrt(3)-1,1/2],
> 

             [-(sqrt(3)+1)/8,1/2,-1/2,-(sqrt(3)+1)/8,1/2,-1/2,-(sqrt(3)+1)/8,-1/2,1/2,3*(sqrt(3)+1)/8]],

> 

            'shape'='symmetric'):
> 

MatrixDecomposition(m, 'method' = 'LDLt'): # this does not work 

Error, (in Student:-NumericalAnalysis:-MatrixDecomposition) a pivot element 0 is encountered, and the entries below it are not all 0; the factorization cannot continue

  
> 

MatrixDecomposition(m([1, 4, 7, 10, 2, 5, 8, 3, 6, 9] $ 2), 'method' = 'LDLt'): # yet this works 
> 

MatrixDecomposition(m([2, 3, 5, 6, 8, 9, 1, 4, 7, 10] $ 2), 'method' = 'LDLt'): # this also works 
> 

randomize(5):
> 

k := 0:
to 1e3 do
        try
                MatrixDecomposition(m(combinat:-randperm(10) $ 2), 'method' = 'LDLt')
        catch :
                k++
        end
od:
k/1e3;

.469

 (1)


 

Download LDL_factorization_robustness.mw

The last instance above suggests that it only works on about half of the inputs (that are equivalent to each other). Although I tried changing the value of Digits, the failure rate remained high. Is Student:-NumericalAnalysis:-MatrixDecomposition not robust enough? 
Please Wait...