put tests in runtests.jl

test/Manifest.toml

+# This file is machine-generated - editing it directly is not advised
+
+[[Base64]]
+uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"
+
+[[Distributed]]
+deps = ["Random", "Serialization", "Sockets"]
+uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b"
+
+[[InteractiveUtils]]
+deps = ["Markdown"]
+uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240"
+
+[[Libdl]]
+uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb"
+
+[[LinearAlgebra]]
+deps = ["Libdl"]
+uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+
+[[Logging]]
+uuid = "56ddb016-857b-54e1-b83d-db4d58db5568"
+
+[[Markdown]]
+deps = ["Base64"]
+uuid = "d6f4376e-aef5-505a-96c1-9c027394607a"
+
+[[Random]]
+deps = ["Serialization"]
+uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+
+[[Serialization]]
+uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b"
+
+[[Sockets]]
+uuid = "6462fe0b-24de-5631-8697-dd941f90decc"
+
+[[Test]]
+deps = ["Distributed", "InteractiveUtils", "Logging", "Random"]
+uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

test/Project.toml

+[deps]
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

test/runtests.jl

+using LinearAlgebra
+using PolynomialMatrixEquations
+using Random
+using Test
+
+undeterminatcase = false
+unstablecas = false
+numberundeterminate = 0
+numberunstable = 0
+
+qz_criterium = 1 + 1e-6
+
+Random.seed!(123)
+ncases = 20
+n=10
+
+for i = 1:ncases
+    println("Test $i")
+    d_orig = randn(n, n)
+    e_orig = randn(n, n)
+    F = schur(e_orig, d_orig)
+    eigenvalues = F.α ./ F.β
+    nstable = count(abs.(eigenvalues) .< 1+1e-6)
+
+    d = copy(d_orig)
+    e = copy(e_orig)
+    ws1 = GsSolverWs(d, e, nstable)
+    gs_solver!(ws1, d, e, nstable, qz_criterium)
+    @test d_orig*[I(nstable); ws1.g2]*ws1.g1 ≈ e_orig*[I(nstable); ws1.g2]
+
+    a0 = Matrix([-e[:, 1:nstable] zeros(n, n-nstable)])
+    a1 = Matrix([d[:, 1:nstable] -e[:, (nstable+1):n]])
+    a2 = Matrix([zeros(n, nstable) d[:, (nstable+1):n]])
+
+    x = zeros(n, n)
+    ws2 = CyclicReductionWs(n)
+    cyclic_reduction!(x, a0, a1, a2, ws2, 1e-8, 50)
+    @test isapprox(a0 + a1*x + a2*x*x, zeros(n, n); atol = 1e-12)
+
+    nstable1 = nstable + 1
+    @test_throws UnstableSystemException gs_solver!(ws1, d, e, nstable + 1, qz_criterium)
+    
+    a0 = Matrix([-e[:, 1:nstable1] zeros(n, n-nstable1)])
+    a1 = Matrix([d[:, 1:nstable1] -e[:, (nstable1+1):n]])
+    a2 = Matrix([zeros(n, nstable1) d[:, (nstable1+1):n]])
+
+    x = zeros(n, n)
+    ws2 = CyclicReductionWs(n)
+
+    @test_throws UnstableSystemException cyclic_reduction!(x, a0, a1, a2, ws2, 1e-8, 50)
+
+    nstable1 = nstable - 1
+    @test_throws UndeterminateSystemException gs_solver!(ws1, d, e, nstable1, qz_criterium)
+
+    a0 = Matrix([-e[:, 1:nstable1] zeros(n, n-nstable1)])
+    a1 = Matrix([d[:, 1:nstable1] -e[:, (nstable1+1):n]])
+    a2 = Matrix([zeros(n, nstable1) d[:, (nstable1+1):n]])
+
+    x = zeros(n, n)
+    ws2 = CyclicReductionWs(n)
+    @test_throws UndeterminateSystemException cyclic_reduction!(x, a0, a1, a2, ws2, 1e-8, 50)
+end

test/testCycleReduction.jl

-using PolynomialMatrixEquations
-using LinearAlgebra
-using Test
-
-n = 3
-ws = CyclicReductionWs(n)
-a0 = [0.5 0 0; 0 0.5 0; 0 0 0];
-a1 = I(n) .+ zeros(n, 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)
-@test ws.info == 0
-@test a0 + a1*x + a2*x*x ≈ zeros(n, n)
-display(x)
-
-cyclic_reduction_check(x,a0,a1,a2,1e-8)
-a0 = [0.5 0 0; 0 1.1 0; 0 0 0];
-a1 .= I(n) .+ zeros(n, n)
-a2 = [0 0 0; 0 0 0; 0 0 0.8]
-cyclic_reduction!(x,a0,a1,a2,ws,1e-8,50)
-@test ws.info == 1
-
-a0 = [0.5 0 0; 0 0.5 0; 0 0 0];
-a1 = I(n) .+ zeros(n, n)
-a2 = [0 0 0; 0 0 0; 0 0 1.2]
-cyclic_reduction!(x,a0,a1,a2,ws,1e-8,50)
-@test ws.info == 2

test/testGeneralizedSchurDecompositionSolver.jl

-using PolynomialMatrixEquations
-using LinearAlgebra
-using Random
-using Test
-
-Random.seed!(132)
-n = 3
-d_org = randn(n,n)
-e_org = randn(n,n)
-
-n1 = count(abs.(eigen(e, d).values) .< 1+1e-6)
-
-ws = GsSolverWs(d, e, n1)
-
-qz_criterium = 1.0 + 1e-6
-d = copy(d_org)
-e = copy(e_org)
-gs_solver!(ws, d , e, n1 , qz_criterium)
-#@test ws.info == 0
-@test d_org*[I(n1); ws.g2]*ws.g1 ≈ e_org*[I(n1); ws.g2]
-
-#=
-cyclic_reduction_check(x,a0,a1,a2,1e-8)
-a0 = [0.5 0 0; 0 1.1 0; 0 0 0];
-a1 .= I(n) .+ zeros(n, n)
-a2 = [0 0 0; 0 0 0; 0 0 0.8]
-cyclic_reduction!(x,a0,a1,a2,ws,1e-8,50)
-@test ws.info == 1
-
-a0 = [0.5 0 0; 0 0.5 0; 0 0 0];
-a1 = I(n) .+ zeros(n, n)
-a2 = [0 0 0; 0 0 0; 0 0 1.2]
-cyclic_reduction!(x,a0,a1,a2,ws,1e-8,50)
-@test ws.info == 2
-=#