Basis for the intersection of vector spaces

How can I get the same functionality in Mathematica?

There exists an old thread called "Intersection of two vector spaces": Intersection of two vector spaces However, I am searching for a general solution which works with any number of vector spaces (like the linalg::intBasis does), not just two.

Moreover, the solution given in this old thread misbehaves in many cases (that I tried).

For three vector spaces I tried getIntersection[V1, getIntersection[V2, V3]] but it returned an incorrect result (where V1, V2 and V3 come from the "Example 1" from the linalg::intBasis web link given above).

Also, if I try getIntersection[V1, V1] it happily returns the {0, 0} vector among the basis vectors.

Additionally, the first comment (by Daniel Lichtblau) in the first "Answer" in this old thread says that one can also use coeffs=ker[[All,1;;n]]; ker.coeffs but the ker.coeffs part often generates an error: Dot::dotsh: Tensors {...} and {...} have incompatible shapes.

(Note: I'm a newcomer here and apparently I have no way to post "comments" in threads which do not belong to me as I get an error saying "You must have 50 reputation to comment" and my "answers" there get deleted. So, I am unable to report problems which I find in them.)

Update (2019.08.07): I have been using the code provided here for several years now. I have found that in some cases, especially when dealing with inexact numerical matrices, it could misbehave. I finally decided to do something about it and, with the help provided in this thread, I came out with this improved function, which may be useful for a casual trespasser who reads this. Note that the Chop function uses a default tolerance of 10.^(-10) but, you can easily change it in the getIntersectionBasis below, e.g.: use Chop[FullSimplify[#], 1.*^-14] (for a tolerance of 10.^(-14)).

getIntersectionBasis[] := {}
getIntersectionBasis[{}] := {}
getIntersectionBasis[{}, __] := {}
getIntersectionBasis[__, {}] := {}
getIntersectionBasis[l1_] := getIntersectionBasis[l1, l1]
getIntersectionBasis[l1_, l2_, l3__] :=
getIntersectionBasis[l1, getIntersectionBasis[l2, l3]]
getIntersectionBasis[l1_, l2_] :=
Catch[With[{ker = FullSimplify[NullSpace[FullSimplify[Transpose[Join[l1, l2]]], ZeroTest->(PossibleZeroQ[Chop[FullSimplify[#]]]&)]]},
If[ker === {}, Throw[{}],
DeleteCases[FullSimplify[RowReduce[FullSimplify[ker[[All, 1 ;; Length[l1]]].l1], ZeroTest->(PossibleZeroQ[Chop[FullSimplify[#]]]&)]], {__?(PossibleZeroQ[Chop[FullSimplify[#]]]&)}]]]
]
• Hi, welcome to Mathematica.SE, please consider taking the tour so you learn the basics of the site. Once you gain enough reputation by making good questions you will be able to vote up and down both questions and answers. When you see good ones, please vote them up by clicking the grey triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. As you receive help, try to give it too, by answering questions in your area of expertise. Oct 26, 2014 at 20:22
• @Wile E. I see that you have posted a question but that it was closed. Please edit your question and explain why it is not a duplicate of the linked one. After you do so I will be glad to reopen it. Oct 27, 2014 at 9:40
• Regarding the inexact case, one might use the Tolerance option in NullSpace. Aug 7, 2019 at 14:29
• Thanks. I did try Tolerance and it never helped (in cases in which I had problems with getIntersectionBasis). Adding ZeroTest solved the problem (at least in the majority of cases). I still have matrices of "exact" complex values for which NullSpace is unable to deliver the result (in these cases I need to apply N and calculate it numerically). Aug 7, 2019 at 15:37

You'd need to row reduce the intersection set to remove linear dependencies it might have. The code below should handle this.

getIntersectionBasis[] := {}
getIntersectionBasis[{}] := {}
getIntersectionBasis[{}, __] := {}
getIntersectionBasis[__, {}] := {}
getIntersectionBasis[l1_] := getIntersectionBasis[l1, l1]
getIntersectionBasis[l1_, l2_, l3__] :=
getIntersectionBasis[getIntersectionBasis[l1, l2], l3]
getIntersectionBasis[l1_, l2_] :=
Catch[With[{ker = NullSpace[Transpose[Join[l1, l2]]]},
If[ker === {}, Throw[{}],
DeleteCases[RowReduce[ker[[All, 1 ;; Length[l1]]].l1], {0 ..}]]]]

Now fill in your favorite example (which actually belongs in the question).

• Thanks for your help. However, it seems to me that, when playing with V1, V2, V3 mentioned above, no matter what I try, the new getIntersectionBasis always returns {{1, 0}, {0, 1}} Oct 27, 2014 at 18:08
• What's wrong with that result? Oct 27, 2014 at 18:40
• Another problem is that when I try the "Example 2" from the linalg::intBasis web link given above, getIntersectionBasis[S1, S2] returns Dot::dotsh: "Tensors {} and {{1,0,1,0},{0,1,0,1}} have incompatible shapes.". Oct 27, 2014 at 19:47
• I'll have a look at that second one. As for the first, well, a basis is not unique. Once we bring in row echelon form (via RowReduce) one is likely to get a result of the sort you are seeing. Oct 27, 2014 at 20:52
• Also, a general rule of thunb: if you want examples tested, include them, in cut-and-paste form, in the query. Not as a link to code posted in a different language. Oct 27, 2014 at 21:08