diff --git a/docs/src/api/solver_errors.md b/docs/src/api/solver_errors.md index 7ce61cdb..49804cdd 100644 --- a/docs/src/api/solver_errors.md +++ b/docs/src/api/solver_errors.md @@ -16,6 +16,9 @@ FullResidual FullResidualRecipe +CompressedResidual + +CompressedResidualRecipe ``` ## Exported Functions diff --git a/src/RLinearAlgebra.jl b/src/RLinearAlgebra.jl index a0c8b9a0..53406a81 100644 --- a/src/RLinearAlgebra.jl +++ b/src/RLinearAlgebra.jl @@ -43,7 +43,7 @@ export LQSolver, LQSolverRecipe # Export SolverError types and functions export SolverError, SolverErrorRecipe export complete_error, compute_error -export FullResidual, FullResidualRecipe +export FullResidual, FullResidualRecipe, CompressedResidual, CompressedResidualRecipe # Export ApproximatorError types and functions export ApproximatorError, ApproximatorErrorRecipe diff --git a/src/Solvers/ErrorMethods.jl b/src/Solvers/ErrorMethods.jl index 93506d21..ee41b877 100644 --- a/src/Solvers/ErrorMethods.jl +++ b/src/Solvers/ErrorMethods.jl @@ -116,3 +116,4 @@ end # Include error method files include("ErrorMethods/full_residual.jl") +include("ErrorMethods/compressed_residual.jl") diff --git a/src/Solvers/ErrorMethods/compressed_residual.jl b/src/Solvers/ErrorMethods/compressed_residual.jl new file mode 100644 index 00000000..7b2cc489 --- /dev/null +++ b/src/Solvers/ErrorMethods/compressed_residual.jl @@ -0,0 +1,54 @@ +""" + CompressedResidual <: ErrorRecipe + +A structure for the compressed residual, `Sb-SAx`. + +# Fields +- None +""" +struct CompressedResidual <: SolverError + +end + +""" + CompressedResidualRecipe <: ErrorRecipe + +A structure for the compressed residual, `Sb-SAx`. + +# Fields +- `residual::AbstractVector`, a container for the compressed residual, `Sb-SAx`. +- `residual_view::SubArray`, a view of the residual container to handle varying compression + sizes. +""" +mutable struct CompressedResidualRecipe{V<:AbstractVector, S<:SubArray} <: SolverErrorRecipe + residual::V + residual_view::S +end + +function complete_error( + error::CompressedResidual, + solver::Solver, + A::AbstractMatrix, + b::AbstractVector +) + residual = zeros(eltype(A), size(b,1)) + residual_view = view(residual, 1:1) + return CompressedResidualRecipe{typeof(residual),typeof(residual_view)}( + residual, + residual_view + ) + +end + +function compute_error( + error::CompressedResidualRecipe, + solver::SolverRecipe, + A::AbstractMatrix, + b::AbstractVector + )::Float64 + rows_s = size(solver.compressor, 1) + error.residual_view = view(error.residual, 1:rows_s) + copyto!(error.residual_view, solver.vec_view) + mul!(error.residual_view, solver.mat_view, solver.solution_vec, -1.0, 1.0) + return norm(error.residual_view) +end diff --git a/test/Solvers/ErrorMethods/compressed_residual.jl b/test/Solvers/ErrorMethods/compressed_residual.jl new file mode 100644 index 00000000..73b34781 --- /dev/null +++ b/test/Solvers/ErrorMethods/compressed_residual.jl @@ -0,0 +1,107 @@ +module compressed_residual_error +using Test, RLinearAlgebra, Random +import LinearAlgebra: mul!, norm +using ..FieldTest +using ..ApproxTol +Random.seed!(1232) + +############################### +# Test Solver Structures +############################### +mutable struct TestSolver <: Solver end + +mutable struct TestSolverRecipe <: SolverRecipe + compressor::AbstractMatrix + vec_view::SubArray + mat_view::SubArray + solution_vec::AbstractVector +end + + +@testset "Compressed Residual" begin + @testset "Compressed Residual: SolverError" begin + # Verify Supertype + @test supertype(CompressedResidual) == SolverError + + # Verify fieldnames and types + @test fieldnames(CompressedResidual) == () + @test fieldtypes(CompressedResidual) == () + # Verify the internal constructor + + end + + @testset "Compressed Residual: SolverErrorRecipe" begin + # Verify Supertype + @test supertype(CompressedResidualRecipe) == SolverErrorRecipe + + # Verify fieldnames and types + @test fieldnames(CompressedResidualRecipe) == (:residual, :residual_view) + @test fieldtypes(CompressedResidualRecipe) == (AbstractVector, SubArray) + end + + @testset "Compressed Residual: Complete error" begin + for type in [Float32, Float64, ComplexF32, ComplexF64] + let n_rows = 4, + n_cols = 3, + comp_dim = 2, + S = ones(type, n_rows, n_cols), + A = ones(type, n_rows, n_cols), + b = ones(type, n_rows), + x = ones(type, n_cols), + error_rec = complete_error(CompressedResidual(), TestSolver(), A, b) + + # Test the type + @test typeof(error_rec) == CompressedResidualRecipe{ + Vector{type}, + SubArray{type, 1, Vector{type}, Tuple{UnitRange{Int64}}, true} + } + # Test type of residual vector + @test eltype(error_rec.residual) == type + # Test residual vector to be all zeros + @test error_rec.residual == zeros(type, n_rows) + end + + end + + end + + @testset "Compressed Residual: Compute Error" begin + for type in [Float32, Float64, ComplexF32, ComplexF64] + let n_rows = 4, + n_cols = 3, + comp_dim = 2, + A = ones(type, n_rows, n_cols), + S = ones(type, comp_dim, n_rows), + b = ones(type, n_rows), + x = ones(type, n_cols), + Sb = S * b, + SA = S * A + solver_rec = TestSolverRecipe( + S, + view(Sb, 1:comp_dim), + view(SA, 1:comp_dim, 1:n_cols), + x + ) + + error_rec = complete_error( + CompressedResidual(), + TestSolver(), + A, + b + ) + + # compute the error value + err_val = compute_error(error_rec, solver_rec, A, b) + # compute the residual + res = Sb - SA * x + # compute norm squared of residual + @test norm(res) ≈ err_val + end + + end + + end + +end + +end