diff --git a/docs/src/reference.md b/docs/src/reference.md
index 0896e1639..5ab99d753 100644
--- a/docs/src/reference.md
+++ b/docs/src/reference.md
@@ -6,7 +6,7 @@
 ```
 
 ```@docs
-Krylov.FloatOrComplex
+Krylov.RealOrComplex
 Krylov.niterations
 Krylov.Aprod
 Krylov.Atprod
diff --git a/src/bicgstab.jl b/src/bicgstab.jl
index 25abd01e6..84648b158 100644
--- a/src/bicgstab.jl
+++ b/src/bicgstab.jl
@@ -21,7 +21,7 @@ export bicgstab, bicgstab!
                           itmax::Int=0, verbose::Int=0, history::Bool=false,
                           ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the square linear system Ax = b using BICGSTAB.
@@ -58,13 +58,13 @@ and `false` otherwise.
 """
 function bicgstab end
 
-function bicgstab(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function bicgstab(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = BicgstabSolver(A, b)
   bicgstab!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function bicgstab(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function bicgstab(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = BicgstabSolver(A, b)
   bicgstab!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -80,7 +80,7 @@ See [`BicgstabSolver`](@ref) for more details about the `solver`.
 """
 function bicgstab! end
 
-function bicgstab!(solver :: BicgstabSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function bicgstab!(solver :: BicgstabSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   bicgstab!(solver, A, b; kwargs...)
   return solver
@@ -89,7 +89,7 @@ end
 function bicgstab!(solver :: BicgstabSolver{T,FC,S}, A, b :: AbstractVector{FC}; c :: AbstractVector{FC}=b,
                    M=I, N=I, atol :: T=√eps(T), rtol :: T=√eps(T),
                    itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-                   ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                   ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("System must be square")
diff --git a/src/bilq.jl b/src/bilq.jl
index 8dbd46c51..35881efe4 100644
--- a/src/bilq.jl
+++ b/src/bilq.jl
@@ -18,7 +18,7 @@ export bilq, bilq!
                       itmax::Int=0, verbose::Int=0, history::Bool=false,
                       callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the square linear system Ax = b using BiLQ.
@@ -46,13 +46,13 @@ and `false` otherwise.
 """
 function bilq end
 
-function bilq(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function bilq(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = BilqSolver(A, b)
   bilq!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function bilq(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function bilq(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = BilqSolver(A, b)
   bilq!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -68,7 +68,7 @@ See [`BilqSolver`](@ref) for more details about the `solver`.
 """
 function bilq! end
 
-function bilq!(solver :: BilqSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function bilq!(solver :: BilqSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   bilq!(solver, A, b; kwargs...)
   return solver
@@ -77,7 +77,7 @@ end
 function bilq!(solver :: BilqSolver{T,FC,S}, A, b :: AbstractVector{FC}; c :: AbstractVector{FC}=b,
                atol :: T=√eps(T), rtol :: T=√eps(T), transfer_to_bicg :: Bool=true,
                itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-               callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("System must be square")
diff --git a/src/bilqr.jl b/src/bilqr.jl
index 479e01319..d50beeff8 100644
--- a/src/bilqr.jl
+++ b/src/bilqr.jl
@@ -18,7 +18,7 @@ export bilqr, bilqr!
                           itmax::Int=0, verbose::Int=0, history::Bool=false,
                           callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Combine BiLQ and QMR to solve adjoint systems.
@@ -48,13 +48,13 @@ and `false` otherwise.
 """
 function bilqr end
 
-function bilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function bilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = BilqrSolver(A, b)
   bilqr!(solver, A, b, c, x0, y0; kwargs...)
   return (solver.x, solver.y, solver.stats)
 end
 
-function bilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function bilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = BilqrSolver(A, b)
   bilqr!(solver, A, b, c; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -71,7 +71,7 @@ See [`BilqrSolver`](@ref) for more details about the `solver`.
 function bilqr! end
 
 function bilqr!(solver :: BilqrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC},
-                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0, y0)
   bilqr!(solver, A, b, c; kwargs...)
   return solver
@@ -80,7 +80,7 @@ end
 function bilqr!(solver :: BilqrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC};
                 atol :: T=√eps(T), rtol :: T=√eps(T), transfer_to_bicg :: Bool=true,
                 itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-                callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("Systems must be square")
diff --git a/src/cg.jl b/src/cg.jl
index 68a6e415d..c2ff02e4f 100644
--- a/src/cg.jl
+++ b/src/cg.jl
@@ -23,7 +23,7 @@ export cg, cg!
                     verbose::Int=0, history::Bool=false,
                     ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 The conjugate gradient method to solve the symmetric linear system Ax = b.
@@ -52,13 +52,13 @@ and `false` otherwise.
 """
 function cg end
 
-function cg(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function cg(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = CgSolver(A, b)
   cg!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function cg(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function cg(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CgSolver(A, b)
   cg!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -74,7 +74,7 @@ See [`CgSolver`](@ref) for more details about the `solver`.
 """
 function cg! end
 
-function cg!(solver :: CgSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function cg!(solver :: CgSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   cg!(solver, A, b; kwargs...)
   return solver
@@ -84,7 +84,7 @@ function cg!(solver :: CgSolver{T,FC,S}, A, b :: AbstractVector{FC};
              M=I, atol :: T=√eps(T), rtol :: T=√eps(T),
              itmax :: Int=0, radius :: T=zero(T), linesearch :: Bool=false,
              verbose :: Int=0, history :: Bool=false,
-             ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+             ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   linesearch && (radius > 0) && error("`linesearch` set to `true` but trust-region radius > 0")
 
diff --git a/src/cg_lanczos.jl b/src/cg_lanczos.jl
index 4e503f09a..f50fdeb8f 100644
--- a/src/cg_lanczos.jl
+++ b/src/cg_lanczos.jl
@@ -19,7 +19,7 @@ export cg_lanczos, cg_lanczos!
                             check_curvature::Bool=false, verbose::Int=0, history::Bool=false,
                             ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 The Lanczos version of the conjugate gradient method to solve the
@@ -48,13 +48,13 @@ and `false` otherwise.
 """
 function cg_lanczos end
 
-function cg_lanczos(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function cg_lanczos(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = CgLanczosSolver(A, b)
   cg_lanczos!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function cg_lanczos(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function cg_lanczos(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CgLanczosSolver(A, b)
   cg_lanczos!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -70,7 +70,7 @@ See [`CgLanczosSolver`](@ref) for more details about the `solver`.
 """
 function cg_lanczos! end
 
-function cg_lanczos!(solver :: CgLanczosSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function cg_lanczos!(solver :: CgLanczosSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   cg_lanczos!(solver, A, b; kwargs...)
   return solver
@@ -79,7 +79,7 @@ end
 function cg_lanczos!(solver :: CgLanczosSolver{T,FC,S}, A, b :: AbstractVector{FC};
                      M=I, atol :: T=√eps(T), rtol :: T=√eps(T), itmax :: Int=0,
                      check_curvature :: Bool=false, verbose :: Int=0, history :: Bool=false,
-                     ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                     ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("System must be square")
diff --git a/src/cg_lanczos_shift.jl b/src/cg_lanczos_shift.jl
index ff873e5b4..2538da90d 100644
--- a/src/cg_lanczos_shift.jl
+++ b/src/cg_lanczos_shift.jl
@@ -21,7 +21,7 @@ export cg_lanczos_shift, cg_lanczos_shift!
                                   verbose::Int=0, history::Bool=false,
                                   ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 The Lanczos version of the conjugate gradient method to solve a family
@@ -39,7 +39,7 @@ and `false` otherwise.
 """
 function cg_lanczos_shift end
 
-function cg_lanczos_shift(A, b :: AbstractVector{FC}, shifts :: AbstractVector{T}; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}}
+function cg_lanczos_shift(A, b :: AbstractVector{FC}, shifts :: AbstractVector{T}; kwargs...) where {T <: Real, FC <: RealOrComplex{T}}
   nshifts = length(shifts)
   solver = CgLanczosShiftSolver(A, b, nshifts)
   cg_lanczos_shift!(solver, A, b, shifts; kwargs...)
@@ -59,7 +59,7 @@ function cg_lanczos_shift!(solver :: CgLanczosShiftSolver{T,FC,S}, A, b :: Abstr
                            M=I, atol :: T=√eps(T), rtol :: T=√eps(T),
                            itmax :: Int=0, check_curvature :: Bool=false,
                            verbose :: Int=0, history :: Bool=false,
-                           ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                           ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("System must be square")
diff --git a/src/cgls.jl b/src/cgls.jl
index a8d6e3c94..e171f1f42 100644
--- a/src/cgls.jl
+++ b/src/cgls.jl
@@ -35,7 +35,7 @@ export cgls, cgls!
                       radius::T=zero(T), itmax::Int=0, verbose::Int=0, history::Bool=false,
                       ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the regularized linear least-squares problem
@@ -63,7 +63,7 @@ and `false` otherwise.
 """
 function cgls end
 
-function cgls(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function cgls(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CglsSolver(A, b)
   cgls!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -81,7 +81,7 @@ function cgls! end
 function cgls!(solver :: CglsSolver{T,FC,S}, A, b :: AbstractVector{FC};
                M=I, λ :: T=zero(T), atol :: T=√eps(T), rtol :: T=√eps(T),
                radius :: T=zero(T), itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/cgne.jl b/src/cgne.jl
index 68039d2de..03e0a430b 100644
--- a/src/cgne.jl
+++ b/src/cgne.jl
@@ -35,7 +35,7 @@ export cgne, cgne!
                       itmax::Int=0, verbose::Int=0, history::Bool=false,
                       ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the consistent linear system
@@ -72,7 +72,7 @@ and `false` otherwise.
 """
 function cgne end
 
-function cgne(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function cgne(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CgneSolver(A, b)
   cgne!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -90,7 +90,7 @@ function cgne! end
 function cgne!(solver :: CgneSolver{T,FC,S}, A, b :: AbstractVector{FC};
                N=I, λ :: T=zero(T), atol :: T=√eps(T), rtol :: T=√eps(T),
                itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/cgs.jl b/src/cgs.jl
index 37a1c4137..e55ea1286 100644
--- a/src/cgs.jl
+++ b/src/cgs.jl
@@ -16,7 +16,7 @@ export cgs, cgs!
                      itmax::Int=0, verbose::Int=0, history::Bool=false,
                      ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the consistent linear system Ax = b using conjugate gradient squared algorithm.
@@ -55,13 +55,13 @@ and `false` otherwise.
 """
 function cgs end
 
-function cgs(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function cgs(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = CgsSolver(A, b)
   cgs!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function cgs(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function cgs(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CgsSolver(A, b)
   cgs!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -77,7 +77,7 @@ See [`CgsSolver`](@ref) for more details about the `solver`.
 """
 function cgs! end
 
-function cgs!(solver :: CgsSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function cgs!(solver :: CgsSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   cgs!(solver, A, b; kwargs...)
   return solver
@@ -86,7 +86,7 @@ end
 function cgs!(solver :: CgsSolver{T,FC,S}, A, b :: AbstractVector{FC}; c :: AbstractVector{FC}=b,
               M=I, N=I, atol :: T=√eps(T), rtol :: T=√eps(T),
               itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-              ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+              ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   m == n || error("System must be square")
diff --git a/src/cr.jl b/src/cr.jl
index c18501b09..d18113ba5 100644
--- a/src/cr.jl
+++ b/src/cr.jl
@@ -20,7 +20,7 @@ export cr, cr!
                     radius::T=zero(T), verbose::Int=0, linesearch::Bool=false, history::Bool=false,
                     ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 A truncated version of Stiefel’s Conjugate Residual method to solve the symmetric linear system Ax = b or the least-squares problem min ‖b - Ax‖.
@@ -51,13 +51,13 @@ and `false` otherwise.
 """
 function cr end
 
-function cr(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function cr(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = CrSolver(A, b)
   cr!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function cr(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function cr(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CrSolver(A, b)
   cr!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -73,7 +73,7 @@ See [`CrSolver`](@ref) for more details about the `solver`.
 """
 function cr! end
 
-function cr!(solver :: CrSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function cr!(solver :: CrSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   cr!(solver, A, b; kwargs...)
   return solver
@@ -82,7 +82,7 @@ end
 function cr!(solver :: CrSolver{T,FC,S}, A, b :: AbstractVector{FC};
              M=I, atol :: T=√eps(T), rtol :: T=√eps(T), γ :: T=√eps(T), itmax :: Int=0,
              radius :: T=zero(T), verbose :: Int=0, linesearch :: Bool=false, history :: Bool=false,
-             ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+             ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   linesearch && (radius > 0) && error("'linesearch' set to 'true' but radius > 0")
   n, m = size(A)
diff --git a/src/craig.jl b/src/craig.jl
index 5759e31df..cbab207b7 100644
--- a/src/craig.jl
+++ b/src/craig.jl
@@ -40,7 +40,7 @@ export craig, craig!
                           verbose::Int=0, transfer_to_lsqr::Bool=false, history::Bool=false,
                           ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Find the least-norm solution of the consistent linear system
@@ -96,7 +96,7 @@ and `false` otherwise.
 """
 function craig end
 
-function craig(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function craig(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CraigSolver(A, b)
   craig!(solver, A, b; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -115,7 +115,7 @@ function craig!(solver :: CraigSolver{T,FC,S}, A, b :: AbstractVector{FC};
                 M=I, N=I, sqd :: Bool=false, λ :: T=zero(T), atol :: T=√eps(T),
                 btol :: T=√eps(T), rtol :: T=√eps(T), conlim :: T=1/√eps(T), itmax :: Int=0,
                 verbose :: Int=0, transfer_to_lsqr :: Bool=false, history :: Bool=false,
-                ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/craigmr.jl b/src/craigmr.jl
index 854e3df98..b124e148d 100644
--- a/src/craigmr.jl
+++ b/src/craigmr.jl
@@ -33,7 +33,7 @@ export craigmr, craigmr!
                             rtol::T=√eps(T), itmax::Int=0, verbose::Int=0, history::Bool=false,
                             ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the consistent linear system
@@ -92,7 +92,7 @@ and `false` otherwise.
 """
 function craigmr end
 
-function craigmr(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function craigmr(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CraigmrSolver(A, b)
   craigmr!(solver, A, b; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -110,7 +110,7 @@ function craigmr! end
 function craigmr!(solver :: CraigmrSolver{T,FC,S}, A, b :: AbstractVector{FC};
                   M=I, N=I, sqd :: Bool=false, λ :: T=zero(T), atol :: T=√eps(T),
                   rtol :: T=√eps(T), itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-                  ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                  ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/crls.jl b/src/crls.jl
index 329b5a5fe..31151dec2 100644
--- a/src/crls.jl
+++ b/src/crls.jl
@@ -27,7 +27,7 @@ export crls, crls!
                       radius::T=zero(T), itmax::Int=0, verbose::Int=0, history::Bool=false,
                       ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the linear least-squares problem
@@ -54,7 +54,7 @@ and `false` otherwise.
 """
 function crls end
 
-function crls(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function crls(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CrlsSolver(A, b)
   crls!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -72,7 +72,7 @@ function crls! end
 function crls!(solver :: CrlsSolver{T,FC,S}, A, b :: AbstractVector{FC};
                M=I, λ :: T=zero(T), atol :: T=√eps(T), rtol :: T=√eps(T),
                radius :: T=zero(T), itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/crmr.jl b/src/crmr.jl
index 3fff12b08..31dace69e 100644
--- a/src/crmr.jl
+++ b/src/crmr.jl
@@ -33,7 +33,7 @@ export crmr, crmr!
                       rtol::T=√eps(T), itmax::Int=0, verbose::Int=0, history::Bool=false,
                       ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the consistent linear system
@@ -70,7 +70,7 @@ and `false` otherwise.
 """
 function crmr end
 
-function crmr(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function crmr(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = CrmrSolver(A, b)
   crmr!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -88,7 +88,7 @@ function crmr! end
 function crmr!(solver :: CrmrSolver{T,FC,S}, A, b :: AbstractVector{FC};
                N=I, λ :: T=zero(T), atol :: T=√eps(T),
                rtol :: T=√eps(T), itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/diom.jl b/src/diom.jl
index 77e73c414..106f59a89 100644
--- a/src/diom.jl
+++ b/src/diom.jl
@@ -17,7 +17,7 @@ export diom, diom!
                       verbose::Int=0, history::Bool=false,
                       ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the consistent linear system Ax = b using DIOM.
@@ -49,13 +49,13 @@ and `false` otherwise.
 """
 function diom end
 
-function diom(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function diom(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = DiomSolver(A, b, memory)
   diom!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function diom(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function diom(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = DiomSolver(A, b, memory)
   diom!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -74,7 +74,7 @@ See [`DiomSolver`](@ref) for more details about the `solver`.
 """
 function diom! end
 
-function diom!(solver :: DiomSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function diom!(solver :: DiomSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   diom!(solver, A, b; kwargs...)
   return solver
@@ -84,7 +84,7 @@ function diom!(solver :: DiomSolver{T,FC,S}, A, b :: AbstractVector{FC};
                M=I, N=I, atol :: T=√eps(T), rtol :: T=√eps(T),
                reorthogonalization :: Bool=false, itmax :: Int=0,
                verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   m == n || error("System must be square")
diff --git a/src/dqgmres.jl b/src/dqgmres.jl
index aa6e245ba..bca2e14cd 100644
--- a/src/dqgmres.jl
+++ b/src/dqgmres.jl
@@ -17,7 +17,7 @@ export dqgmres, dqgmres!
                          verbose::Int=0, history::Bool=false,
                          ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the consistent linear system Ax = b using DQGMRES.
@@ -49,13 +49,13 @@ and `false` otherwise.
 """
 function dqgmres end
 
-function dqgmres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function dqgmres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = DqgmresSolver(A, b, memory)
   dqgmres!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function dqgmres(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function dqgmres(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = DqgmresSolver(A, b, memory)
   dqgmres!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -74,7 +74,7 @@ See [`DqgmresSolver`](@ref) for more details about the `solver`.
 """
 function dqgmres! end
 
-function dqgmres!(solver :: DqgmresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function dqgmres!(solver :: DqgmresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   dqgmres!(solver, A, b; kwargs...)
   return solver
@@ -84,7 +84,7 @@ function dqgmres!(solver :: DqgmresSolver{T,FC,S}, A, b :: AbstractVector{FC};
                   M=I, N=I, atol :: T=√eps(T), rtol :: T=√eps(T),
                   reorthogonalization :: Bool=false, itmax :: Int=0,
                   verbose :: Int=0, history :: Bool=false,
-                  ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                  ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   m == n || error("System must be square")
diff --git a/src/fgmres.jl b/src/fgmres.jl
index 635e7241e..fd34396c8 100644
--- a/src/fgmres.jl
+++ b/src/fgmres.jl
@@ -17,7 +17,7 @@ export fgmres, fgmres!
                         restart::Bool=false, verbose::Int=0, history::Bool=false,
                         ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the linear system Ax = b using FGMRES.
@@ -52,13 +52,13 @@ and `false` otherwise.
 """
 function fgmres end
 
-function fgmres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function fgmres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = FgmresSolver(A, b, memory)
   fgmres!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function fgmres(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function fgmres(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = FgmresSolver(A, b, memory)
   fgmres!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -77,7 +77,7 @@ See [`FgmresSolver`](@ref) for more details about the `solver`.
 """
 function fgmres! end
 
-function fgmres!(solver :: FgmresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function fgmres!(solver :: FgmresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   fgmres!(solver, A, b; kwargs...)
   return solver
@@ -87,7 +87,7 @@ function fgmres!(solver :: FgmresSolver{T,FC,S}, A, b :: AbstractVector{FC};
                 M=I, N=I, atol :: T=√eps(T), rtol :: T=√eps(T),
                 reorthogonalization :: Bool=false, itmax :: Int=0,
                 restart :: Bool=false, verbose :: Int=0, history :: Bool=false,
-                ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   m == n || error("System must be square")
diff --git a/src/fom.jl b/src/fom.jl
index 95bcc97d1..da6e31ae8 100644
--- a/src/fom.jl
+++ b/src/fom.jl
@@ -17,7 +17,7 @@ export fom, fom!
                      restart::Bool=false, verbose::Int=0, history::Bool=false,
                      ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the linear system Ax = b using FOM.
@@ -46,13 +46,13 @@ and `false` otherwise.
 """
 function fom end
 
-function fom(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function fom(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = FomSolver(A, b, memory)
   fom!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function fom(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function fom(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = FomSolver(A, b, memory)
   fom!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -71,7 +71,7 @@ See [`FomSolver`](@ref) for more details about the `solver`.
 """
 function fom! end
 
-function fom!(solver :: FomSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function fom!(solver :: FomSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   fom!(solver, A, b; kwargs...)
   return solver
@@ -81,7 +81,7 @@ function fom!(solver :: FomSolver{T,FC,S}, A, b :: AbstractVector{FC};
               M=I, N=I, atol :: T=√eps(T), rtol :: T=√eps(T),
               reorthogonalization :: Bool=false, itmax :: Int=0,
               restart :: Bool=false, verbose :: Int=0, history :: Bool=false,
-              ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+              ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   m == n || error("System must be square")
diff --git a/src/gmres.jl b/src/gmres.jl
index b145b512b..beb659312 100644
--- a/src/gmres.jl
+++ b/src/gmres.jl
@@ -17,7 +17,7 @@ export gmres, gmres!
                        restart::Bool=false, verbose::Int=0, history::Bool=false,
                        ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the linear system Ax = b using GMRES.
@@ -46,13 +46,13 @@ and `false` otherwise.
 """
 function gmres end
 
-function gmres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function gmres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = GmresSolver(A, b, memory)
   gmres!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function gmres(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function gmres(A, b :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = GmresSolver(A, b, memory)
   gmres!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -71,7 +71,7 @@ See [`GmresSolver`](@ref) for more details about the `solver`.
 """
 function gmres! end
 
-function gmres!(solver :: GmresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function gmres!(solver :: GmresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   gmres!(solver, A, b; kwargs...)
   return solver
@@ -81,7 +81,7 @@ function gmres!(solver :: GmresSolver{T,FC,S}, A, b :: AbstractVector{FC};
                 M=I, N=I, atol :: T=√eps(T), rtol :: T=√eps(T),
                 reorthogonalization :: Bool=false, itmax :: Int=0,
                 restart :: Bool=false, verbose :: Int=0, history :: Bool=false,
-                ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   m == n || error("System must be square")
diff --git a/src/gpmr.jl b/src/gpmr.jl
index 528bd522d..8efcb830c 100644
--- a/src/gpmr.jl
+++ b/src/gpmr.jl
@@ -19,7 +19,7 @@ export gpmr, gpmr!
                          verbose::Int=0, history::Bool=false,
                          ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 GPMR solves the unsymmetric partitioned linear system
@@ -73,13 +73,13 @@ and `false` otherwise.
 """
 function gpmr end
 
-function gpmr(A, B, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function gpmr(A, B, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = GpmrSolver(A, b, memory)
   gpmr!(solver, A, B, b, c, x0, y0; kwargs...)
   return (solver.x, solver.y, solver.stats)
 end
 
-function gpmr(A, B, b :: AbstractVector{FC}, c :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: FloatOrComplex
+function gpmr(A, B, b :: AbstractVector{FC}, c :: AbstractVector{FC}; memory :: Int=20, kwargs...) where FC <: RealOrComplex
   solver = GpmrSolver(A, b, memory)
   gpmr!(solver, A, B, b, c; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -99,7 +99,7 @@ See [`GpmrSolver`](@ref) for more details about the `solver`.
 function gpmr! end
 
 function gpmr!(solver :: GpmrSolver{T,FC,S}, A, B, b :: AbstractVector{FC}, c :: AbstractVector{FC},
-                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0, y0)
   gpmr!(solver, A, B, b, c; kwargs...)
   return solver
@@ -110,7 +110,7 @@ function gpmr!(solver :: GpmrSolver{T,FC,S}, A, B, b :: AbstractVector{FC}, c ::
                gsp :: Bool=false, reorthogonalization :: Bool=false,
                itmax :: Int=0, λ :: FC=one(FC), μ :: FC=one(FC),
                verbose :: Int=0, history::Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   s, t = size(B)
diff --git a/src/krylov_processes.jl b/src/krylov_processes.jl
index 49a30f4b3..4b55a85c6 100644
--- a/src/krylov_processes.jl
+++ b/src/krylov_processes.jl
@@ -18,7 +18,7 @@ export hermitian_lanczos, nonhermitian_lanczos, arnoldi, golub_kahan, saunders_s
 
 * C. Lanczos, [*An Iteration Method for the Solution of the Eigenvalue Problem of Linear Differential and Integral Operators*](https://doi.org/10.6028/jres.045.026), Journal of Research of the National Bureau of Standards, 45(4), pp. 225--280, 1950.
 """
-function hermitian_lanczos(A, b::AbstractVector{FC}, k::Int) where FC <: FloatOrComplex
+function hermitian_lanczos(A, b::AbstractVector{FC}, k::Int) where FC <: RealOrComplex
   n, m = size(A)
   R = real(FC)
   S = ktypeof(b)
@@ -89,7 +89,7 @@ end
 
 * C. Lanczos, [*An Iteration Method for the Solution of the Eigenvalue Problem of Linear Differential and Integral Operators*](https://doi.org/10.6028/jres.045.026), Journal of Research of the National Bureau of Standards, 45(4), pp. 225--280, 1950.
 """
-function nonhermitian_lanczos(A, b::AbstractVector{FC}, c::AbstractVector{FC}, k::Int) where FC <: FloatOrComplex
+function nonhermitian_lanczos(A, b::AbstractVector{FC}, c::AbstractVector{FC}, k::Int) where FC <: RealOrComplex
   n, m = size(A)
   Aᴴ = A'
   S = ktypeof(b)
@@ -178,7 +178,7 @@ end
 
 * W. E. Arnoldi, [*The principle of minimized iterations in the solution of the matrix eigenvalue problem*](https://doi.org/10.1090/qam/42792), Quarterly of Applied Mathematics, 9, pp. 17--29, 1951.
 """
-function arnoldi(A, b::AbstractVector{FC}, k::Int) where FC <: FloatOrComplex
+function arnoldi(A, b::AbstractVector{FC}, k::Int) where FC <: RealOrComplex
   n, m = size(A)
   S = ktypeof(b)
   M = vector_to_matrix(S)
@@ -238,7 +238,7 @@ end
 
 * G. H. Golub and W. Kahan, [*Calculating the Singular Values and Pseudo-Inverse of a Matrix*](https://doi.org/10.1137/0702016), SIAM Journal on Numerical Analysis, 2(2), pp. 225--224, 1965.
 """
-function golub_kahan(A, b::AbstractVector{FC}, k::Int) where FC <: FloatOrComplex
+function golub_kahan(A, b::AbstractVector{FC}, k::Int) where FC <: RealOrComplex
   n, m = size(A)
   R = real(FC)
   Aᴴ = A'
@@ -313,7 +313,7 @@ end
 
 * M. A. Saunders, H. D. Simon, and E. L. Yip, [*Two Conjugate-Gradient-Type Methods for Unsymmetric Linear Equations*](https://doi.org/10.1137/0725052), SIAM Journal on Numerical Analysis, 25(4), pp. 927--940, 1988.
 """
-function saunders_simon_yip(A, b::AbstractVector{FC}, c::AbstractVector{FC}, k::Int) where FC <: FloatOrComplex
+function saunders_simon_yip(A, b::AbstractVector{FC}, c::AbstractVector{FC}, k::Int) where FC <: RealOrComplex
   n, m = size(A)
   Aᴴ = A'
   S = ktypeof(b)
@@ -404,7 +404,7 @@ end
 
 * A. Montoison and D. Orban, [*GPMR: An Iterative Method for Unsymmetric Partitioned Linear Systems*](https://dx.doi.org/10.13140/RG.2.2.24069.68326), Cahier du GERAD G-2021-62, GERAD, Montréal, 2021.
 """
-function montoison_orban(A, B, b::AbstractVector{FC}, c::AbstractVector{FC}, k::Int) where FC <: FloatOrComplex
+function montoison_orban(A, B, b::AbstractVector{FC}, c::AbstractVector{FC}, k::Int) where FC <: RealOrComplex
   n, m = size(A)
   S = ktypeof(b)
   M = vector_to_matrix(S)
diff --git a/src/krylov_solvers.jl b/src/krylov_solvers.jl
index f94efd2f9..6bb32a45f 100644
--- a/src/krylov_solvers.jl
+++ b/src/krylov_solvers.jl
@@ -1817,7 +1817,7 @@ end
 
 Statistics of `solver` are displayed if `show_stats` is set to true.
 """
-function show(io :: IO, solver :: KrylovSolver{T,FC,S}; show_stats :: Bool=true) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function show(io :: IO, solver :: KrylovSolver{T,FC,S}; show_stats :: Bool=true) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   workspace = typeof(solver)
   name_solver = workspace.name.wrapper
   l1 = max(length(string(name_solver)), 10)  # length("warm_start") = 10
diff --git a/src/krylov_stats.jl b/src/krylov_stats.jl
index a662fa0a0..0f0cb3b4b 100644
--- a/src/krylov_stats.jl
+++ b/src/krylov_stats.jl
@@ -225,7 +225,7 @@ for f in ["Simple", "Lsmr", "Lanczos", "LanczosShift", "Symmlq", "Adjoint", "LNL
       end
       s *=  " " * field_name * ":"
       statfield = getfield(stats, field)
-      if isa(statfield, AbstractVector) && eltype(statfield) <: Union{Missing, AbstractFloat}
+      if isa(statfield, AbstractVector) && eltype(statfield) <: Union{Missing, Real}
         s *= @sprintf " %s\n" vec2str(statfield)
       else
         s *= @sprintf " %s\n" statfield
diff --git a/src/krylov_utils.jl b/src/krylov_utils.jl
index b0a180dd3..61f8b779d 100644
--- a/src/krylov_utils.jl
+++ b/src/krylov_utils.jl
@@ -1,8 +1,8 @@
 """
-    FloatOrComplex{T}
-Union type of `T` and `Complex{T}` where T is an `AbstractFloat`.
+    RealOrComplex{T}
+Union type of `T` and `Complex{T}` where T is a `Real`.
 """
-const FloatOrComplex{T} = Union{T, Complex{T}} where T <: AbstractFloat
+const RealOrComplex{T} = Union{T, Complex{T}} where T <: Real
 
 """
     (c, s, ρ) = sym_givens(a, b)
@@ -13,7 +13,7 @@ Given `a` and `b` reals, return `(c, s, ρ)` such that
     [ c  s ] [ a ] = [ ρ ]
     [ s -c ] [ b ] = [ 0 ].
 """
-function sym_givens(a :: T, b :: T) where T <: AbstractFloat
+function sym_givens(a :: T, b :: T) where T <: Real
   #
   # Modeled after the corresponding Matlab function by M. A. Saunders and S.-C. Choi.
   # http://www.stanford.edu/group/SOL/dissertations/sou-cheng-choi-thesis.pdf
@@ -57,7 +57,7 @@ c real and (s, ρ) complexes such that
     [ c   s ] [ a ] = [ ρ ]
     [ s̅  -c ] [ b ] = [ 0 ].
 """
-function sym_givens(a :: Complex{T}, b :: Complex{T}) where T <: AbstractFloat
+function sym_givens(a :: Complex{T}, b :: Complex{T}) where T <: Real
   #
   # Modeled after the corresponding Fortran function by M. A. Saunders and S.-C. Choi.
   # A. Montoison, Montreal, March 2020.
@@ -92,8 +92,8 @@ function sym_givens(a :: Complex{T}, b :: Complex{T}) where T <: AbstractFloat
   return (c, s, ρ)
 end
 
-sym_givens(a :: Complex{T}, b :: T) where T <: AbstractFloat = sym_givens(a, Complex{T}(b))
-sym_givens(a :: T, b :: Complex{T}) where T <: AbstractFloat = sym_givens(Complex{T}(a), b)
+sym_givens(a :: Complex{T}, b :: T) where T <: Real = sym_givens(a, Complex{T}(b))
+sym_givens(a :: T, b :: Complex{T}) where T <: Real = sym_givens(Complex{T}(a), b)
 
 """
     roots = roots_quadratic(q₂, q₁, q₀; nitref)
@@ -107,7 +107,7 @@ cancellation. Optionally, `nitref` steps of iterative refinement
 may be performed to improve accuracy. By default, `nitref=1`.
 """
 function roots_quadratic(q₂ :: T, q₁ :: T, q₀ :: T;
-                         nitref :: Int=1) where T <: AbstractFloat
+                         nitref :: Int=1) where T <: Real
   # Case where q(x) is linear.
   if q₂ == zero(T)
     if q₁ == zero(T)
@@ -159,7 +159,7 @@ Display an array in the form
 
 with (ndisp - 1)/2 elements on each side.
 """
-function vec2str(x :: AbstractVector{T}; ndisp :: Int=7) where T <: Union{AbstractFloat, Missing}
+function vec2str(x :: AbstractVector{T}; ndisp :: Int=7) where T <: Union{Real, Missing}
   n = length(x)
   if n ≤ ndisp
     ndisp = n
@@ -262,30 +262,30 @@ mulorldiv!(y, P, x, ldiv::Bool) = ldiv ? ldiv!(y, P, x) : mul!(y, P, x)
 
 kdot(n :: Integer, x :: Vector{T}, dx :: Integer, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasReal = BLAS.dot(n, x, dx, y, dy)
 kdot(n :: Integer, x :: Vector{T}, dx :: Integer, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasComplex = BLAS.dotc(n, x, dx, y, dy)
-kdot(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: FloatOrComplex = dot(x, y)
+kdot(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: RealOrComplex = dot(x, y)
 
-kdotr(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: AbstractFloat = kdot(n, x, dx, y, dy)
-kdotr(n :: Integer, x :: AbstractVector{Complex{T}}, dx :: Integer, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = real(kdot(n, x, dx, y, dy))
+kdotr(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: Real = kdot(n, x, dx, y, dy)
+kdotr(n :: Integer, x :: AbstractVector{Complex{T}}, dx :: Integer, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: Real = real(kdot(n, x, dx, y, dy))
 
 knrm2(n :: Integer, x :: Vector{T}, dx :: Integer) where T <: BLAS.BlasFloat = BLAS.nrm2(n, x, dx)
-knrm2(n :: Integer, x :: AbstractVector{T}, dx :: Integer) where T <: FloatOrComplex = norm(x)
+knrm2(n :: Integer, x :: AbstractVector{T}, dx :: Integer) where T <: RealOrComplex = norm(x)
 
 kscal!(n :: Integer, s :: T, x :: Vector{T}, dx :: Integer) where T <: BLAS.BlasFloat = BLAS.scal!(n, s, x, dx)
-kscal!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer) where T <: FloatOrComplex = (x .*= s)
-kscal!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer) where T <: AbstractFloat = kscal!(n, Complex{T}(s), x, dx)
+kscal!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer) where T <: RealOrComplex = (x .*= s)
+kscal!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer) where T <: Real = kscal!(n, Complex{T}(s), x, dx)
 
 kaxpy!(n :: Integer, s :: T, x :: Vector{T}, dx :: Integer, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasFloat = BLAS.axpy!(n, s, x, dx, y, dy)
-kaxpy!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: FloatOrComplex = axpy!(s, x, y)
-kaxpy!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = kaxpy!(n, Complex{T}(s), x, dx, y, dy)
+kaxpy!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: RealOrComplex = axpy!(s, x, y)
+kaxpy!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: Real = kaxpy!(n, Complex{T}(s), x, dx, y, dy)
 
 kaxpby!(n :: Integer, s :: T, x :: Vector{T}, dx :: Integer, t :: T, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasFloat = BLAS.axpby!(n, s, x, dx, t, y, dy)
-kaxpby!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer, t :: T, y :: AbstractVector{T}, dy :: Integer) where T <: FloatOrComplex = axpby!(s, x, t, y)
-kaxpby!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: Complex{T}, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = kaxpby!(n, Complex{T}(s), x, dx, t, y, dy)
-kaxpby!(n :: Integer, s :: Complex{T}, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: T, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = kaxpby!(n, s, x, dx, Complex{T}(t), y, dy)
-kaxpby!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: T, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = kaxpby!(n, Complex{T}(s), x, dx, Complex{T}(t), y, dy)
+kaxpby!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer, t :: T, y :: AbstractVector{T}, dy :: Integer) where T <: RealOrComplex = axpby!(s, x, t, y)
+kaxpby!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: Complex{T}, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: Real = kaxpby!(n, Complex{T}(s), x, dx, t, y, dy)
+kaxpby!(n :: Integer, s :: Complex{T}, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: T, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: Real = kaxpby!(n, s, x, dx, Complex{T}(t), y, dy)
+kaxpby!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: T, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: Real = kaxpby!(n, Complex{T}(s), x, dx, Complex{T}(t), y, dy)
 
 kcopy!(n :: Integer, x :: Vector{T}, dx :: Integer, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasFloat = BLAS.blascopy!(n, x, dx, y, dy)
-kcopy!(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: FloatOrComplex = copyto!(y, x)
+kcopy!(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: RealOrComplex = copyto!(y, x)
 
 # the macros are just for readability, so we don't have to write the increments (always equal to 1)
 macro kdot(n, x, y)
@@ -344,7 +344,7 @@ If `flip` is set to `true`, `σ1` and `σ2` are computed such that
 
     ‖x - σi d‖ = radius, i = 1, 2.
 """
-function to_boundary(n :: Int, x :: AbstractVector{T}, d :: AbstractVector{T}, radius :: T; flip :: Bool=false, xNorm2 :: T=zero(T), dNorm2 :: T=zero(T)) where T <: FloatOrComplex
+function to_boundary(n :: Int, x :: AbstractVector{T}, d :: AbstractVector{T}, radius :: T; flip :: Bool=false, xNorm2 :: T=zero(T), dNorm2 :: T=zero(T)) where T <: RealOrComplex
   radius > 0 || error("radius must be positive")
 
   # ‖d‖² σ² + (xᴴd + dᴴx) σ + (‖x‖² - Δ²).
diff --git a/src/lnlq.jl b/src/lnlq.jl
index db0a7c951..5d41fbeeb 100644
--- a/src/lnlq.jl
+++ b/src/lnlq.jl
@@ -31,7 +31,7 @@ export lnlq, lnlq!
                          transfer_to_craig::Bool=true, verbose::Int=0, history::Bool=false,
                          ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Find the least-norm solution of the consistent linear system
@@ -88,7 +88,7 @@ and `false` otherwise.
 """
 function lnlq end
 
-function lnlq(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function lnlq(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = LnlqSolver(A, b)
   lnlq!(solver, A, b; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -107,7 +107,7 @@ function lnlq!(solver :: LnlqSolver{T,FC,S}, A, b :: AbstractVector{FC};
                M=I, N=I, sqd :: Bool=false, λ :: T=zero(T), σ :: T=zero(T),
                atol :: T=√eps(T), rtol :: T=√eps(T), etolx :: T=√eps(T), etoly :: T=√eps(T), itmax :: Int=0,
                transfer_to_craig :: Bool=true, verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/lslq.jl b/src/lslq.jl
index d43d4a089..f640d75bb 100644
--- a/src/lslq.jl
+++ b/src/lslq.jl
@@ -31,7 +31,7 @@ export lslq, lslq!
                       conlim::T=1/√eps(T), verbose::Int=0, history::Bool=false,
                       ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the regularized linear least-squares problem
@@ -137,7 +137,7 @@ and `false` otherwise.
 """
 function lslq end
 
-function lslq(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: FloatOrComplex
+function lslq(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: RealOrComplex
   solver = LslqSolver(A, b, window=window)
   lslq!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -158,7 +158,7 @@ function lslq!(solver :: LslqSolver{T,FC,S}, A, b :: AbstractVector{FC};
                utol :: T=√eps(T), itmax :: Int=0, σ :: T=zero(T),
                transfer_to_lsqr :: Bool=false, conlim :: T=1/√eps(T),
                verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/lsmr.jl b/src/lsmr.jl
index 79d2543fb..41b847f55 100644
--- a/src/lsmr.jl
+++ b/src/lsmr.jl
@@ -36,7 +36,7 @@ export lsmr, lsmr!
                       history::Bool=false, ldiv::Bool=false,
                       callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the regularized linear least-squares problem
@@ -94,7 +94,7 @@ and `false` otherwise.
 """
 function lsmr end
 
-function lsmr(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: FloatOrComplex
+function lsmr(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: RealOrComplex
   solver = LsmrSolver(A, b, window=window)
   lsmr!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -115,7 +115,7 @@ function lsmr!(solver :: LsmrSolver{T,FC,S}, A, b :: AbstractVector{FC};
                atol :: T=zero(T), rtol :: T=zero(T),
                etol :: T=√eps(T), itmax :: Int=0, conlim :: T=1/√eps(T),
                radius :: T=zero(T), verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/lsqr.jl b/src/lsqr.jl
index e4973bd38..d02609ba7 100644
--- a/src/lsqr.jl
+++ b/src/lsqr.jl
@@ -35,7 +35,7 @@ export lsqr, lsqr!
                       radius::T=zero(T), verbose::Int=0, history::Bool=false,
                       ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the regularized linear least-squares problem
@@ -89,7 +89,7 @@ and `false` otherwise.
 """
 function lsqr end
 
-function lsqr(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: FloatOrComplex
+function lsqr(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: RealOrComplex
   solver = LsqrSolver(A, b, window=window)
   lsqr!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -110,7 +110,7 @@ function lsqr!(solver :: LsqrSolver{T,FC,S}, A, b :: AbstractVector{FC};
                atol :: T=zero(T), rtol :: T=zero(T),
                etol :: T=√eps(T), itmax :: Int=0, conlim :: T=1/√eps(T),
                radius :: T=zero(T), verbose :: Int=0, history :: Bool=false,
-               ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+               ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/minres.jl b/src/minres.jl
index c95048bbc..4b5a224bc 100644
--- a/src/minres.jl
+++ b/src/minres.jl
@@ -32,7 +32,7 @@ export minres, minres!
                         history::Bool=false, ldiv::Bool=false,
                         callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the shifted linear least-squares problem
@@ -70,13 +70,13 @@ and `false` otherwise.
 """
 function minres end
 
-function minres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; window :: Int=5, kwargs...) where FC <: FloatOrComplex
+function minres(A, b :: AbstractVector{FC}, x0 :: AbstractVector; window :: Int=5, kwargs...) where FC <: RealOrComplex
   solver = MinresSolver(A, b, window=window)
   minres!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function minres(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: FloatOrComplex
+function minres(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: RealOrComplex
   solver = MinresSolver(A, b, window=window)
   minres!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -92,7 +92,7 @@ See [`MinresSolver`](@ref) for more details about the `solver`.
 """
 function minres! end
 
-function minres!(solver :: MinresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function minres!(solver :: MinresSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   minres!(solver, A, b; kwargs...)
   return solver
@@ -102,7 +102,7 @@ function minres!(solver :: MinresSolver{T,FC,S}, A, b :: AbstractVector{FC};
                  M=I, λ :: T=zero(T), atol :: T=√eps(T)/100, rtol :: T=√eps(T)/100, 
                  ratol :: T=zero(T), rrtol :: T=zero(T), etol :: T=√eps(T),
                  itmax :: Int=0, conlim :: T=1/√eps(T), verbose :: Int=0,
-                 history :: Bool=false, ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                 history :: Bool=false, ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("System must be square")
diff --git a/src/minres_qlp.jl b/src/minres_qlp.jl
index cb70754b8..d069e3e4e 100644
--- a/src/minres_qlp.jl
+++ b/src/minres_qlp.jl
@@ -23,7 +23,7 @@ export minres_qlp, minres_qlp!
                             verbose::Int=0, history::Bool=false,
                             ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 MINRES-QLP is the only method based on the Lanczos process that returns the minimum-norm
@@ -51,13 +51,13 @@ and `false` otherwise.
 """
 function minres_qlp end
 
-function minres_qlp(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function minres_qlp(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = MinresQlpSolver(A, b)
   minres_qlp!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function minres_qlp(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function minres_qlp(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = MinresQlpSolver(A, b)
   minres_qlp!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -73,7 +73,7 @@ See [`MinresQlpSolver`](@ref) for more details about the `solver`.
 """
 function minres_qlp! end
 
-function minres_qlp!(solver :: MinresQlpSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function minres_qlp!(solver :: MinresQlpSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   minres_qlp!(solver, A, b; kwargs...)
   return solver
@@ -83,7 +83,7 @@ function minres_qlp!(solver :: MinresQlpSolver{T,FC,S}, A, b :: AbstractVector{F
                      M=I, atol :: T=√eps(T), rtol :: T=√eps(T),
                      ctol :: T=√eps(T), λ ::T=zero(T), itmax :: Int=0,
                      verbose :: Int=0, history :: Bool=false,
-                     ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                     ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("System must be square")
diff --git a/src/qmr.jl b/src/qmr.jl
index fe0fab65c..f41888206 100644
--- a/src/qmr.jl
+++ b/src/qmr.jl
@@ -26,7 +26,7 @@ export qmr, qmr!
                      itmax::Int=0, verbose::Int=0, history::Bool=false,
                      callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the square linear system Ax = b using QMR.
@@ -52,13 +52,13 @@ and `false` otherwise.
 """
 function qmr end
 
-function qmr(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function qmr(A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = QmrSolver(A, b)
   qmr!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function qmr(A, b :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function qmr(A, b :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = QmrSolver(A, b)
   qmr!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -74,7 +74,7 @@ See [`QmrSolver`](@ref) for more details about the `solver`.
 """
 function qmr! end
 
-function qmr!(solver :: QmrSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function qmr!(solver :: QmrSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   qmr!(solver, A, b; kwargs...)
   return solver
@@ -83,7 +83,7 @@ end
 function qmr!(solver :: QmrSolver{T,FC,S}, A, b :: AbstractVector{FC}; c :: AbstractVector{FC}=b,
               atol :: T=√eps(T), rtol :: T=√eps(T),
               itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-              callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+              callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   n, m = size(A)
   m == n || error("System must be square")
diff --git a/src/symmlq.jl b/src/symmlq.jl
index efbd751aa..e1013599c 100644
--- a/src/symmlq.jl
+++ b/src/symmlq.jl
@@ -20,7 +20,7 @@ export symmlq, symmlq!
                         verbose::Int=0, history::Bool=false,
                         ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the shifted linear system
@@ -50,13 +50,13 @@ and `false` otherwise.
 """
 function symmlq end
 
-function symmlq(A, b :: AbstractVector{FC}, x0 :: AbstractVector; window :: Int=5, kwargs...) where FC <: FloatOrComplex
+function symmlq(A, b :: AbstractVector{FC}, x0 :: AbstractVector; window :: Int=5, kwargs...) where FC <: RealOrComplex
   solver = SymmlqSolver(A, b, window=window)
   symmlq!(solver, A, b, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function symmlq(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: FloatOrComplex
+function symmlq(A, b :: AbstractVector{FC}; window :: Int=5, kwargs...) where FC <: RealOrComplex
   solver = SymmlqSolver(A, b, window=window)
   symmlq!(solver, A, b; kwargs...)
   return (solver.x, solver.stats)
@@ -72,7 +72,7 @@ See [`SymmlqSolver`](@ref) for more details about the `solver`.
 """
 function symmlq! end
 
-function symmlq!(solver :: SymmlqSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+function symmlq!(solver :: SymmlqSolver{T,FC,S}, A, b :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   symmlq!(solver, A, b; kwargs...)
   return solver
@@ -83,7 +83,7 @@ function symmlq!(solver :: SymmlqSolver{T,FC,S}, A, b :: AbstractVector{FC};
                  λest :: T=zero(T), atol :: T=√eps(T), rtol :: T=√eps(T),
                  etol :: T=√eps(T), itmax :: Int=0, conlim :: T=1/√eps(T),
                  verbose :: Int=0, history :: Bool=false,
-                 ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                 ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   m == n || error("System must be square")
diff --git a/src/tricg.jl b/src/tricg.jl
index 8d0a41ce3..45fb1a9e3 100644
--- a/src/tricg.jl
+++ b/src/tricg.jl
@@ -19,7 +19,7 @@ export tricg, tricg!
                           verbose::Int=0, history::Bool=false,
                           ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 TriCG solves the symmetric linear system
@@ -68,13 +68,13 @@ and `false` otherwise.
 """
 function tricg end
 
-function tricg(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function tricg(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = TricgSolver(A, b)
   tricg!(solver, A, b, c, x0, y0; kwargs...)
   return (solver.x, solver.y, solver.stats)
 end
 
-function tricg(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function tricg(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = TricgSolver(A, b)
   tricg!(solver, A, b, c; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -91,7 +91,7 @@ See [`TricgSolver`](@ref) for more details about the `solver`.
 function tricg! end
 
 function tricg!(solver :: TricgSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC},
-                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0, y0)
   tricg!(solver, A, b, c; kwargs...)
   return solver
@@ -102,7 +102,7 @@ function tricg!(solver :: TricgSolver{T,FC,S}, A, b :: AbstractVector{FC}, c ::
                 spd :: Bool=false, snd :: Bool=false, flip :: Bool=false,
                 τ :: T=one(T), ν :: T=-one(T), itmax :: Int=0,
                 verbose :: Int=0, history :: Bool=false,
-                ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/trilqr.jl b/src/trilqr.jl
index 60663ff55..b2ab489d6 100644
--- a/src/trilqr.jl
+++ b/src/trilqr.jl
@@ -18,7 +18,7 @@ export trilqr, trilqr!
                            itmax::Int=0, verbose::Int=0, history::Bool=false,
                            callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Combine USYMLQ and USYMQR to solve adjoint systems.
@@ -47,13 +47,13 @@ and `false` otherwise.
 """
 function trilqr end
 
-function trilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function trilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = TrilqrSolver(A, b)
   trilqr!(solver, A, b, c, x0, y0; kwargs...)
   return (solver.x, solver.y, solver.stats)
 end
 
-function trilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function trilqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = TrilqrSolver(A, b)
   trilqr!(solver, A, b, c; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -70,7 +70,7 @@ See [`TrilqrSolver`](@ref) for more details about the `solver`.
 function trilqr! end
 
 function trilqr!(solver :: TrilqrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC},
-                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0, y0)
   trilqr!(solver, A, b, c; kwargs...)
   return solver
@@ -79,7 +79,7 @@ end
 function trilqr!(solver :: TrilqrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC};
                  atol :: T=√eps(T), rtol :: T=√eps(T), transfer_to_usymcg :: Bool=true,
                  itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-                 callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                 callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/trimr.jl b/src/trimr.jl
index 041a5ffff..f9cbef608 100644
--- a/src/trimr.jl
+++ b/src/trimr.jl
@@ -19,7 +19,7 @@ export trimr, trimr!
                           verbose::Int=0, history::Bool=false,
                           ldiv::Bool=false, callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 TriMR solves the symmetric linear system
@@ -68,13 +68,13 @@ and `false` otherwise.
 """
 function trimr end
 
-function trimr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function trimr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = TrimrSolver(A, b)
   trimr!(solver, A, b, c, x0, y0; kwargs...)
   return (solver.x, solver.y, solver.stats)
 end
 
-function trimr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function trimr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = TrimrSolver(A, b)
   trimr!(solver, A, b, c; kwargs...)
   return (solver.x, solver.y, solver.stats)
@@ -91,7 +91,7 @@ See [`TrimrSolver`](@ref) for more details about the `solver`.
 function trimr! end
 
 function trimr!(solver :: TrimrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC},
-                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                x0 :: AbstractVector, y0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0, y0)
   trimr!(solver, A, b, c; kwargs...)
   return solver
@@ -102,7 +102,7 @@ function trimr!(solver :: TrimrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c ::
                 spd :: Bool=false, snd :: Bool=false, flip :: Bool=false, sp :: Bool=false,
                 τ :: T=one(T), ν :: T=-one(T), itmax :: Int=0,
                 verbose :: Int=0, history :: Bool=false,
-                ldiv :: Bool=false, callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                ldiv :: Bool=false, callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/usymlq.jl b/src/usymlq.jl
index acec8d77e..a5178d892 100644
--- a/src/usymlq.jl
+++ b/src/usymlq.jl
@@ -25,7 +25,7 @@ export usymlq, usymlq!
                         itmax::Int=0, verbose::Int=0, history::Bool=false,
                         callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the linear system Ax = b using the USYMLQ method.
@@ -58,13 +58,13 @@ and `false` otherwise.
 """
 function usymlq end
 
-function usymlq(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function usymlq(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = UsymlqSolver(A, b)
   usymlq!(solver, A, b, c, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function usymlq(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function usymlq(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = UsymlqSolver(A, b)
   usymlq!(solver, A, b, c; kwargs...)
   return (solver.x, solver.stats)
@@ -81,7 +81,7 @@ See [`UsymlqSolver`](@ref) for more details about the `solver`.
 function usymlq! end
 
 function usymlq!(solver :: UsymlqSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC},
-                 x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                 x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   usymlq!(solver, A, b, c; kwargs...)
   return solver
@@ -90,7 +90,7 @@ end
 function usymlq!(solver :: UsymlqSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC};
                  atol :: T=√eps(T), rtol :: T=√eps(T), transfer_to_usymcg :: Bool=true,
                  itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-                 callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                 callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/src/usymqr.jl b/src/usymqr.jl
index 13c19efa8..71843e165 100644
--- a/src/usymqr.jl
+++ b/src/usymqr.jl
@@ -25,7 +25,7 @@ export usymqr, usymqr!
                         itmax::Int=0, verbose::Int=0, history::Bool=false,
                         callback=solver->false)
 
-`T` is an `AbstractFloat` such as `Float32`, `Float64` or `BigFloat`.
+`T` is a `Real` such as `Float32`, `Float64` or `BigFloat`.
 `FC` is `T` or `Complex{T}`.
 
 Solve the linear system Ax = b using USYMQR.
@@ -55,13 +55,13 @@ and `false` otherwise.
 """
 function usymqr end
 
-function usymqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: FloatOrComplex
+function usymqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}, x0 :: AbstractVector; kwargs...) where FC <: RealOrComplex
   solver = UsymqrSolver(A, b)
   usymqr!(solver, A, b, c, x0; kwargs...)
   return (solver.x, solver.stats)
 end
 
-function usymqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: FloatOrComplex
+function usymqr(A, b :: AbstractVector{FC}, c :: AbstractVector{FC}; kwargs...) where FC <: RealOrComplex
   solver = UsymqrSolver(A, b)
   usymqr!(solver, A, b, c; kwargs...)
   return (solver.x, solver.stats)
@@ -78,7 +78,7 @@ See [`UsymqrSolver`](@ref) for more details about the `solver`.
 function usymqr! end
 
 function usymqr!(solver :: UsymqrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC},
-                 x0 :: AbstractVector; kwargs...) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                 x0 :: AbstractVector; kwargs...) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
   warm_start!(solver, x0)
   usymqr!(solver, A, b, c; kwargs...)
   return solver
@@ -87,7 +87,7 @@ end
 function usymqr!(solver :: UsymqrSolver{T,FC,S}, A, b :: AbstractVector{FC}, c :: AbstractVector{FC};
                  atol :: T=√eps(T), rtol :: T=√eps(T),
                  itmax :: Int=0, verbose :: Int=0, history :: Bool=false,
-                 callback = solver -> false) where {T <: AbstractFloat, FC <: FloatOrComplex{T}, S <: DenseVector{FC}}
+                 callback = solver -> false) where {T <: Real, FC <: RealOrComplex{T}, S <: DenseVector{FC}}
 
   m, n = size(A)
   length(b) == m || error("Inconsistent problem size")
diff --git a/test/gpu/intel.jl b/test/gpu/intel.jl
index e65a2bf47..3ef78b7d0 100644
--- a/test/gpu/intel.jl
+++ b/test/gpu/intel.jl
@@ -4,7 +4,7 @@ include("gpu.jl")
 
 import Krylov.kdot
 # https://github.com/JuliaGPU/GPUArrays.jl/pull/427
-function kdot(n :: Integer, x :: oneVector{T}, dx :: Integer, y :: oneVector{T}, dy :: Integer) where T <: Krylov.FloatOrComplex
+function kdot(n :: Integer, x :: oneVector{T}, dx :: Integer, y :: oneVector{T}, dy :: Integer) where T <: Krylov.RealOrComplex
   return mapreduce(dot, +, x, y)
 end
 
diff --git a/test/gpu/metal.jl b/test/gpu/metal.jl
index 9fd24392a..904c681e4 100644
--- a/test/gpu/metal.jl
+++ b/test/gpu/metal.jl
@@ -8,7 +8,7 @@ const MtlMatrix{T} = MtlArray{T,2}
 
 # https://github.com/JuliaGPU/GPUArrays.jl/pull/427
 import Krylov.kdot
-function kdot(n :: Integer, x :: MtlVector{T}, dx :: Integer, y :: MtlVector{T}, dy :: Integer) where T <: Krylov.FloatOrComplex
+function kdot(n :: Integer, x :: MtlVector{T}, dx :: Integer, y :: MtlVector{T}, dy :: Integer) where T <: Krylov.RealOrComplex
   return mapreduce(dot, +, x, y)
 end
 
diff --git a/test/test_stats.jl b/test/test_stats.jl
index 4289a78a3..54c241a76 100644
--- a/test/test_stats.jl
+++ b/test/test_stats.jl
@@ -15,7 +15,7 @@
   @test strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
   check_reset(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 
   stats = Krylov.LsmrStats(0, true, true, Float64[1.0], Float64[2.0], Float64(3.0), Float64(4.0), Float64(5.0), Float64(6.0), Float64(7.0), "t")
   io = IOBuffer()
@@ -37,7 +37,7 @@
   @test strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
   check_reset(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 
   stats = Krylov.LanczosStats(0, true, Float64[3.0], true, NaN, NaN, "t")
   io = IOBuffer()
@@ -55,7 +55,7 @@
   @test strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
   check_reset(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 
   stats = Krylov.LanczosShiftStats(0, true, [Float64[0.9, 0.5], Float64[0.6, 0.4, 0.1]], BitVector([false, true]), NaN, NaN, "t")
   io = IOBuffer()
@@ -70,9 +70,9 @@
   ‖A‖F: NaN
   κ₂(A): NaN
   status: t"""
-  @test (VERSION < v"1.5") || strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
+  strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 
   stats = Krylov.SymmlqStats(0, true, Float64[4.0], Union{Float64,Missing}[5.0, missing], Float64[6.0], Union{Float64,Missing}[7.0, missing], NaN, NaN, "t")
   io = IOBuffer()
@@ -92,7 +92,7 @@
   @test strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
   check_reset(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 
   stats = Krylov.AdjointStats(0, true, true, Float64[8.0], Float64[9.0], "t")
   io = IOBuffer()
@@ -109,7 +109,7 @@
   @test strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
   check_reset(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 
   stats = Krylov.LNLQStats(0, true, Float64[10.0], false, Float64[11.0], Float64[12.0], "t")
   io = IOBuffer()
@@ -127,7 +127,7 @@
   @test strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
   check_reset(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 
   stats = Krylov.LSLQStats(0, true, false, Float64[13.0], Float64[14.0], Float64[15.0], false, Float64[16.0], Float64[17.0], "t")
   io = IOBuffer()
@@ -148,5 +148,5 @@
   @test strip.(split(chomp(showed), "\n")) == strip.(split(chomp(expected), "\n"))
   Krylov.reset!(stats)
   check_reset(stats)
-  @test (VERSION < v"1.5") || (@allocated Krylov.reset!(stats)) == 0
+  (@allocated Krylov.reset!(stats)) == 0
 end