document and test colleague_matrix

Author: stevengj

README.md

@@ -32,7 +32,8 @@ the latter case, `c(y)` returns a vector of interpolants, and one can use `chebj
 to compute the tuple `(c(y), J(y))` of the interpolant and its Jacobian matrix at `y`.
 
 There is also a function `roots(c)` to return an array of the real roots of `c` on the interval `[lb, ub]`,
-computed efficiently by a ["colleague"-matrix method (combined with recursive domain decomposition)](https://epubs.siam.org/doi/10.1137/1.9781611975949.ch18).
+computed efficiently by a ["colleague"-matrix method (combined with recursive domain decomposition)](https://epubs.siam.org/doi/10.1137/1.9781611975949.ch18).  One can also obtain the colleague matrix (whose eigenvalues are the roots)
+directly via `colleague_matrix(c)` (which can be useful to examine complex roots and roots slightly outside `[lb, ub]`).
 
 ### Regression from arbitrary points
 

test/runtests.jl

@@ -29,6 +29,7 @@ Random.seed!(314159) # make chainrules tests deterministic
         @test roots(chebinterp(x -> x - 0.1, 10, 0,1)) ≈ [0.1]
         @test roots(chebinterp(x -> (x - 0.1)*(x-0.2), 10, 0,1)) ≈ [0.1, 0.2]
         @test roots(chebinterp(cos, 10000, 0, 1000))/pi ≈ (0:317) .+ T(0.5)
+        @test colleague_matrix(chebinterp(x -> (x - 0.1)*(x-0.2), 10, 0,1)) ≈ T[0.5 -0.24; 0.5 -0.2]
     end
 end