Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

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.

share|improve this question
1  
It's just a guess but perhaps it's because all elements in the matrix Flatten[Inverse[IdentityMatrix[3] - #]] aren't unique? In that case SubsetQ returns false since Intersection[{1,1,2},{1,2}] is not {1,1,2} but {1,2}. The length of the list is two, not three. –  Pickett Aug 18 '13 at 23:11
    
@Anon, thanks, you definitely found a bug. However that wasn't my main problem. I edited the above to include the error messages I'm getting, which might be helpful. –  Twiffy Aug 19 '13 at 1:01

2 Answers 2

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.

share|improve this answer
up vote 1 down vote accepted

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]
share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.