How do I create matrices with arbitrary list restrictions?

Question

I want to randomly generate a square matrix of dimension n with entries in the list StartingEntries, and also satisfying that for any such matrix M, the matrix Inverse[IdentityMatrix[n] - M] has entries coming from the list EndingEntries.

Below is my (not-working) code for the case dimension 3. It's the last line that keeps generating error messages. I'd like to know why. (Also if anyone can think of a better way of solving this problem, I'm all ears!)

(*generate 10,000 possible matrices.*)
mat = {{#1, #2, #3}, {#4, #5, #6}, {#7, #8, #9}} &
candidates = Apply[mat, RandomChoice[StartingEntries, {10000, 9}], 1]

(*get all matrices M so that I-M is invertible*)
invertibleCandidates = Select[candidates, Det[IdentityMatrix[3] - #] != 0 &]

So far so good. Now we want to make sure that the inverted matrix has the correct entries. I couldn't find a boolean Subset function in Mathematica, so I made my own:

SubsetQ[x_, y_] := (Length[Intersection[x, y]] == Length[x])

Now test to see if the inverse has entries lying inside EndingEntries. This is the line that keeps failing on me and I don't know why.

newerCan = Select[invertibleCandidates, SubsetQ[Flatten[Inverse[IdentityMatrix[3] - #]], EndingEntries] &]

Edit 1: SubSetQ modified to be

SubsetQ[x_, y_] := (Length[Intersection[DeleteDuplicates[x], DeleteDuplicates[y]]] == Length[DeleteDuplicates[x]])

Edit 2: The error messages I'm getting are

DeleteDuplicates::normal: Nonatomic expression expected at position 1 in DeleteDuplicates[nums].
Intersection::heads: Heads DeleteDuplicates and List at positions 2 and 1 are expected to be the same.
Inverse::luc: "Result for Inverse of badly conditioned matrix {{0.6,-0.7,-0.2},{-0.4,0.6,-0.1},{-0.2,-0.3,1.}} may contain significant numerical errors."
Inverse::sing: "Matrix {{0.7,-0.8,0.},{0.,0.9,-0.9},{-0.7,0.,0.8}} is singular."

I'm assuming that the first three error messages are because I'm not formatting the input to SubSetQ in the appropriate list format, but I don't know how to fix it, and I don't know why I'd be getting the "singular" error message.

Answer 1

First evaluate Inverse[{{1, 1}, {2, 2}}] and you will see that the expression is not evaluated since the matrix is singular. This causes the first two error messages, because you're trying to execute DeleteDuplicates on an unevaluated (atomic) expression.

Because some matrices are badly conditioned, the values you get from those matrices won't be exact and hence they won't match the elements in EndingEntries, or they will match but it will be incorrect. That's one problem but it's not fatal, it just affect the detectability of certain matrices.

The error message Inverse::sing: "Matrix {{0.7,-0.8,0.},{0.,0.9,-0.9},{-0.7,0.,0.8}} is singular." is exactly what you wanted to avoid by checking the determinant. The reason it didn't work can be seen to be that:

Det[{{0.7, -0.8, 0.}, {0., 0.9, -0.9}, {-0.7, 0., 0.8}}]

7.99361*10^-17

That is, because of problems with precision, Mathematica reports this matrix to be non-singular. But in reality:

MatrixRank[{{0.7, -0.8, 0.}, {0., 0.9, -0.9}, {-0.7, 0., 0.8}}]

2

So that its determinant ought to have been zero, because this matrix is evidently singular.

You can either introduce a tolerance level and replace the condition !=0 with Abs[]>tolerance or something like that. Or you can use this test: MatrixRank[m]=n. This should ensure that only truly non-singular/invertible matrices are selected into invertibleCandidates.

Answer 2

Turns out I had two issues: one was that Mathematica was incorrectly identifying my matrices as "badly conditioned" just because the entries were decimals. I changed them to fractions (of the same value) and those error messages went away. Second, there was something wonky about my list formatting, and I still don't understand exactly what. However, the following code is now working perfectly -- it's essentially the same as what's above, except that I keep everything as a 1-dimensional list, and use Partition to convert to a matrix when necessary, instead of keeping everything as a matrix and flattening to a list when necessary.

(*create our list of candidates*)
candidates = RandomChoice[StartingEntries, {10000, 9}]

(*identify the matrices M where I-M is invertible*)
invertible = Select[candidates, Det[IdentityMatrix[3] - Partition[#, 3]] != 0 &]

(*define the SubsetQ and ValidMatrix functions*)
SubsetQ[x_, y_] := Length[Intersection[DeleteDuplicates[x], DeleteDuplicates[y]]] == Length[DeleteDuplicates[x]]
ValidMatrix[x_] := SubsetQ[Flatten[Inverse[IdentityMatrix[3] - Partition[x, 3]]], EndingEntries]

(*select the ones where Inverse[I-M] have the appropriate entries*)
finalChoices = Select[invertible, ValidMatrix]

(*get one randomly*)
Partition[RandomChoice[finalChoices], 3]