Matlab - generate random nonsingular triangular matrix

Question

How can I generate a random triangular matrix? (upper and lower)

Usually I use rand(n) but if I try tril(rand(n)) it will be singular and I don't want that.

Answer 1

your answer is correct:

 A=tril(rand(n))

you can check that this matrix is not singular using

 rcond(A)>eps

or

 min(svd(A))>eps

and verifying that the smallest singular value is larger than eps, or any other numerical tolerance that is relevant to your needs. (the code will return 1 or 0). For n>50 you'll start to approach singular matrices.

Here's a small analysis of how the matrix approaches singularity with it's size...

Answer 2

Ok, here are the thoughts on the triangular matrix singularity. Determinant of the triangular matrix is what determines singularity, because it is going into denominator when inverse matrix is build. The property of the triangular matrix is such that determinant equals to the product of diagonal elements.

So for matrix NxN we on the diagonal have the product of i.i.d. U(0,1) numbers. Obviously, determinant will decrease as N increases due to the fact that all number are <1 and more you have, the less value the product (aka determinant) will have.

It is interesting to check that for det=X₁X₂...*X_N the mean value will go down as 2^-N, because each term in product is U(0,1) with mean of 1/2, and they are all i.i.d. Alternative check would be to compute mean from product PDF (see https://math.stackexchange.com/questions/659254/product-distribution-of-two-uniform-distribution-what-about-3-or-more), and, indeed, it will give you exactly the same result, 2^-N. Variance of the determinant could be computed as well, as second momentum minus mean squared, and it is equal to (3^-N-4^-N).

Please note, those are mean values, which you could expect on average, say, if you sample 10⁶ triangular matrices with N=100, compute their determinant and average it, you should find it to be pretty close to 2^-100.

That's where the problem lies. On average, triangular random matrix goes to singularity exponentially with the grows of N. 2^-10 is about equal to 1/1,000. 2^-20 is about equal to 1/1,000,000. For N=100 it should be, on average, around 10^-30 or so, which makes whole exercise moot.

Unfortunately, I cannot offer anything beyond this simple analysis.

Answer 3

If you want a well-conditioned random triangular matrix you could take the triangular part of A² with A = rand(n). So triu(A * A) for any size n is well conditioned but of course has complexity O(n³) for the matrix-matrix multiplication.