From 31057a56821198a4c3e7976e5e8a459b9155ac04 Mon Sep 17 00:00:00 2001 From: Thore Kockerols Date: Thu, 7 Nov 2024 10:16:44 +0000 Subject: [PATCH] default schur solvers for qme --- src/MacroModelling.jl | 14 +++++++------- src/algorithms/quadratic_matrix_equation.jl | 6 ++++-- src/get_functions.jl | 2 +- src/perturbation.jl | 6 +++--- 4 files changed, 15 insertions(+), 13 deletions(-) diff --git a/src/MacroModelling.jl b/src/MacroModelling.jl index 89197a75..86b953c8 100644 --- a/src/MacroModelling.jl +++ b/src/MacroModelling.jl @@ -3864,7 +3864,7 @@ function calculate_second_order_stochastic_steady_state(parameters::Vector{M}, 𝓂::ℳ; verbose::Bool = false, pruning::Bool = false, - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :doubling, timer::TimerOutput = TimerOutput(), tol::AbstractFloat = 1e-12)::Tuple{Vector{M}, Bool, Vector{M}, M, AbstractMatrix{M}, SparseMatrixCSC{M}, AbstractMatrix{M}, SparseMatrixCSC{M}} where M @@ -4199,7 +4199,7 @@ function calculate_third_order_stochastic_steady_state( parameters::Vector{M}, 𝓂::ℳ; verbose::Bool = false, pruning::Bool = false, - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :bicgstab, timer::TimerOutput = TimerOutput(), tol::AbstractFloat = 1e-12)::Tuple{Vector{M}, Bool, Vector{M}, M, AbstractMatrix{M}, SparseMatrixCSC{M}, SparseMatrixCSC{M}, AbstractMatrix{M}, SparseMatrixCSC{M}, SparseMatrixCSC{M}} where M @@ -7061,7 +7061,7 @@ function get_relevant_steady_state_and_state_update(::Val{:second_order}, parameter_values::Vector{S}, 𝓂::ℳ, tol::AbstractFloat; - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :doubling, verbose::Bool = false, timer::TimerOutput = TimerOutput()) where S <: Real @@ -7092,7 +7092,7 @@ function get_relevant_steady_state_and_state_update(::Val{:pruned_second_order}, parameter_values::Vector{S}, 𝓂::ℳ, tol::AbstractFloat; - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :doubling, verbose::Bool = false, timer::TimerOutput = TimerOutput())::Tuple{timings, Vector{S}, Union{Matrix{S},Vector{AbstractMatrix{S}}}, Vector{Vector{S}}, Bool} where S <: Real @@ -7125,7 +7125,7 @@ function get_relevant_steady_state_and_state_update(::Val{:third_order}, parameter_values::Vector{S}, 𝓂::ℳ, tol::AbstractFloat; - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :bicgstab, verbose::Bool = false, timer::TimerOutput = TimerOutput()) where S <: Real @@ -7156,7 +7156,7 @@ function get_relevant_steady_state_and_state_update(::Val{:pruned_third_order}, parameter_values::Vector{S}, 𝓂::ℳ, tol::AbstractFloat; - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :bicgstab, verbose::Bool = false, timer::TimerOutput = TimerOutput())::Tuple{timings, Vector{S}, Union{Matrix{S},Vector{AbstractMatrix{S}}}, Vector{Vector{S}}, Bool} where S <: Real @@ -7187,7 +7187,7 @@ function get_relevant_steady_state_and_state_update(::Val{:first_order}, parameter_values::Vector{S}, 𝓂::ℳ, tol::AbstractFloat; - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :gmres, verbose::Bool = false, timer::TimerOutput = TimerOutput())::Tuple{timings, Vector{S}, Union{Matrix{S},Vector{AbstractMatrix{S}}}, Vector{Vector{Float64}}, Bool} where S <: Real diff --git a/src/algorithms/quadratic_matrix_equation.jl b/src/algorithms/quadratic_matrix_equation.jl index ee820c94..334ad659 100644 --- a/src/algorithms/quadratic_matrix_equation.jl +++ b/src/algorithms/quadratic_matrix_equation.jl @@ -11,7 +11,7 @@ function solve_quadratic_matrix_equation(A::AbstractMatrix{R}, C::AbstractMatrix{R}, T::timings; initial_guess::AbstractMatrix{R} = zeros(0,0), - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, timer::TimerOutput = TimerOutput(), tol::AbstractFloat = 1e-10, # 1e-14 is too tight verbose::Bool = false) where R <: Real @@ -69,6 +69,8 @@ function solve_quadratic_matrix_equation(A::AbstractMatrix{R}, end end + if (reached_tol > tol) println("QME failed: $reached_tol") end + return sol, reached_tol < tol end @@ -507,7 +509,7 @@ function solve_quadratic_matrix_equation(A::AbstractMatrix{ℱ.Dual{Z,S,N}}, T::timings; initial_guess::AbstractMatrix{<:Real} = zeros(0,0), tol::AbstractFloat = 1e-10, - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, timer::TimerOutput = TimerOutput(), verbose::Bool = false) where {Z,S,N} # unpack: AoS -> SoA diff --git a/src/get_functions.jl b/src/get_functions.jl index b77bac5c..e34ca952 100644 --- a/src/get_functions.jl +++ b/src/get_functions.jl @@ -3015,7 +3015,7 @@ function get_loglikelihood(𝓂::ℳ, presample_periods::Int = 0, initial_covariance::Symbol = :theoretical, filter_algorithm::Symbol = :LagrangeNewton, - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, sylvester_algorithm::Symbol = :bicgstab, tol::AbstractFloat = 1e-12, timer::TimerOutput = TimerOutput(), diff --git a/src/perturbation.jl b/src/perturbation.jl index 78e0f8bc..a231ff46 100644 --- a/src/perturbation.jl +++ b/src/perturbation.jl @@ -1,6 +1,6 @@ function calculate_first_order_solution(∇₁::Matrix{Float64}; T::timings, - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, verbose::Bool = false, initial_guess::AbstractMatrix{<:AbstractFloat} = zeros(0,0), timer::TimerOutput = TimerOutput())::Tuple{Matrix{Float64}, Matrix{Float64}, Bool} @@ -119,7 +119,7 @@ end function rrule(::typeof(calculate_first_order_solution), ∇₁::Matrix{Float64}; T::timings, - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, verbose::Bool = false, initial_guess::AbstractMatrix{<:AbstractFloat} = zeros(0,0), timer::TimerOutput = TimerOutput()) @@ -274,7 +274,7 @@ end function calculate_first_order_solution(∇₁::Matrix{ℱ.Dual{Z,S,N}}; T::timings, - quadratic_matrix_equation_solver::Symbol = :doubling, + quadratic_matrix_equation_solver::Symbol = :schur, verbose::Bool = false, initial_guess::AbstractMatrix{<:AbstractFloat} = zeros(0,0), timer::TimerOutput = TimerOutput())::Tuple{Matrix{ℱ.Dual{Z,S,N}}, Matrix{Float64}, Bool} where {Z,S,N}