WORK IN PROGRESS

This package solves matrix polynomial equations in the form

A0 + A1*X + A2*X*X = 0

or

 -    -       -    -
 | I  |       | I  |
D| -  | X1 = E| -  |
 | X2 |       | X2 |
 -    -       -    -

where matrices X1 and X2 contain columns from the solution matrix X.

Two algorithms are provided

Cyclic reduction

cyclic_reduction!(x::Array{Float64},a0::Array{Float64},a1::Array{Float64},a2::Array{Float64},ws::CyclicReductionWs, cvg_tol::Float64, max_it::Int64)

Solve the quadratic matrix equation a0 + a1x + a2x*x = 0, using the cyclic reduction method from Bini et al. (???).

The solution is returned in matrix x. In case of nonconvergency, x is set to NaN and UndeterminateSystemExcpetion or UnstableSystemException is thrown

Example

using CyclicReduction
n = 3
ws = CyclicReductionWs(n)
a0 = [0.5 0 0; 0 0.5 0; 0 0 0];
a1 = eye(n)
a2 = [0 0 0; 0 0 0; 0 0 0.8]
x = zeros(n,n)


cyclic_reduction!(x,a0,a1,a2,ws,1e-8,50)



## General Schur Decomposition

gs_solver!(ws::GsSolverWs,d::Matrix{Float64},e::Matrix{Float64},n1::Int64,qz_criterium)

The solution is returned in ``ws.g1`` and ``ws.g2``