Skip to content

Commit

Permalink
Test Krylov macros
Browse files Browse the repository at this point in the history
amontoison committed Sep 9, 2022
1 parent 73e95a7 commit fc8677b
Showing 2 changed files with 246 additions and 193 deletions.
151 changes: 75 additions & 76 deletions src/krylov_utils.jl
Original file line number Diff line number Diff line change
@@ -175,6 +175,49 @@ function to_boundary(x :: Vector{T}, d :: Vector{T},
return roots # `σ1` and `σ2`
end

"""
s = vec2str(x; ndisp)
Display an array in the form
[ -3.0e-01 -5.1e-01 1.9e-01 ... -2.3e-01 -4.4e-01 2.4e-01 ]
with (ndisp - 1)/2 elements on each side.
"""
function vec2str(x :: AbstractVector{T}; ndisp :: Int=7) where T <: Union{AbstractFloat, Missing}
n = length(x)
if n ndisp
ndisp = n
nside = n
else
nside = max(1, div(ndisp - 1, 2))
end
s = "["
i = 1
while i nside
if x[i] !== missing
s *= @sprintf("%8.1e ", x[i])
else
s *= " ✗✗✗✗ "
end
i += 1
end
if i div(n, 2)
s *= "... "
end
i = max(i, n - nside + 1)
while i n
if x[i] !== missing
s *= @sprintf("%8.1e ", x[i])
else
s *= " ✗✗✗✗ "
end
i += 1
end
s *= "]"
return s
end

"""
S = ktypeof(v)
@@ -209,76 +252,75 @@ end
Create an AbstractVector of storage type `S` of length `n` only composed of zero.
"""
@inline kzeros(S, n) = fill!(S(undef, n), zero(eltype(S)))
kzeros(S, n) = fill!(S(undef, n), zero(eltype(S)))

"""
v = kones(S, n)
Create an AbstractVector of storage type `S` of length `n` only composed of one.
"""
@inline kones(S, n) = fill!(S(undef, n), one(eltype(S)))
kones(S, n) = fill!(S(undef, n), one(eltype(S)))

@inline allocate_if(bool, solver, v, S, n) = bool && isempty(solver.:($v)) && (solver.:($v) = S(undef, n))
allocate_if(bool, solver, v, S, n) = bool && isempty(solver.:($v)) && (solver.:($v) = S(undef, n))

@inline kdisplay(iter, verbose) = (verbose > 0) && (mod(iter, verbose) == 0)
kdisplay(iter, verbose) = (verbose > 0) && (mod(iter, verbose) == 0)

@inline mulorldiv!(y, P, x, ldiv::Bool) = ldiv ? ldiv!(y, P, x) : mul!(y, P, x)
mulorldiv!(y, P, x, ldiv::Bool) = ldiv ? ldiv!(y, P, x) : mul!(y, P, x)

@inline krylov_dot(n :: Integer, x :: Vector{T}, dx :: Integer, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasReal = BLAS.dot(n, x, dx, y, dy)
@inline krylov_dot(n :: Integer, x :: Vector{T}, dx :: Integer, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasComplex = BLAS.dotc(n, x, dx, y, dy)
@inline krylov_dot(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: Number = dot(x, y)
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)

@inline krylov_dotr(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: AbstractFloat = krylov_dot(n, x, dx, y, dy)
@inline krylov_dotr(n :: Integer, x :: AbstractVector{Complex{T}}, dx :: Integer, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = real(krylov_dot(n, x, dx, y, dy))
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))

@inline krylov_norm2(n :: Integer, x :: Vector{T}, dx :: Integer) where T <: BLAS.BlasFloat = BLAS.nrm2(n, x, dx)
@inline krylov_norm2(n :: Integer, x :: AbstractVector{T}, dx :: Integer) where T <: Number = norm(x)
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)

@inline krylov_scal!(n :: Integer, s :: T, x :: Vector{T}, dx :: Integer) where T <: BLAS.BlasFloat = BLAS.scal!(n, s, x, dx)
@inline krylov_scal!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer) where T <: Number = (x .*= s)
@inline krylov_scal!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer) where T <: AbstractFloat = krylov_scal!(n, Complex{T}(s), x, dx)
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)

@inline krylov_axpy!(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)
@inline krylov_axpy!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: Number = axpy!(s, x, y)
@inline krylov_axpy!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = krylov_axpy!(n, Complex{T}(s), x, dx, y, dy)
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)

@inline krylov_axpby!(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)
@inline krylov_axpby!(n :: Integer, s :: T, x :: AbstractVector{T}, dx :: Integer, t :: T, y :: AbstractVector{T}, dy :: Integer) where T <: Number = axpby!(s, x, t, y)
@inline krylov_axpby!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: Complex{T}, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = krylov_axpby!(n, Complex{T}(s), x, dx, t, y, dy)
@inline krylov_axpby!(n :: Integer, s :: Complex{T}, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: T, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = krylov_axpby!(n, s, x, dx, Complex{T}(t), y, dy)
@inline krylov_axpby!(n :: Integer, s :: T, x :: AbstractVector{Complex{T}}, dx :: Integer, t :: T, y :: AbstractVector{Complex{T}}, dy :: Integer) where T <: AbstractFloat = krylov_axpby!(n, Complex{T}(s), x, dx, Complex{T}(t), 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)

@inline krylov_copy!(n :: Integer, x :: Vector{T}, dx :: Integer, y :: Vector{T}, dy :: Integer) where T <: BLAS.BlasFloat = BLAS.blascopy!(n, x, dx, y, dy)
@inline krylov_copy!(n :: Integer, x :: AbstractVector{T}, dx :: Integer, y :: AbstractVector{T}, dy :: Integer) where T <: Number = copyto!(y, x)
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)

# the macros are just for readability, so we don't have to write the increments (always equal to 1)

macro kdot(n, x, y)
return esc(:(krylov_dot($n, $x, 1, $y, 1)))
return esc(:(Krylov.kdot($n, $x, 1, $y, 1)))
end

macro kdotr(n, x, y)
return esc(:(krylov_dotr($n, $x, 1, $y, 1)))
return esc(:(Krylov.kdotr($n, $x, 1, $y, 1)))
end

macro knrm2(n, x)
return esc(:(krylov_norm2($n, $x, 1)))
return esc(:(Krylov.knrm2($n, $x, 1)))
end

macro kscal!(n, s, x)
return esc(:(krylov_scal!($n, $s, $x, 1)))
return esc(:(Krylov.kscal!($n, $s, $x, 1)))
end

macro kaxpy!(n, s, x, y)
return esc(:(krylov_axpy!($n, $s, $x, 1, $y, 1)))
return esc(:(Krylov.kaxpy!($n, $s, $x, 1, $y, 1)))
end

macro kaxpby!(n, s, x, t, y)
return esc(:(krylov_axpby!($n, $s, $x, 1, $t, $y, 1)))
return esc(:(Krylov.kaxpby!($n, $s, $x, 1, $t, $y, 1)))
end

macro kcopy!(n, x, y)
return esc(:(krylov_copy!($n, $x, 1, $y, 1)))
return esc(:(Krylov.kcopy!($n, $x, 1, $y, 1)))
end

macro kswap(x, y)
@@ -292,46 +334,3 @@ end
macro kref!(n, x, y, c, s)
return esc(:(reflect!($x, $y, $c, $s)))
end

"""
s = vec2str(x; ndisp)
Display an array in the form
[ -3.0e-01 -5.1e-01 1.9e-01 ... -2.3e-01 -4.4e-01 2.4e-01 ]
with (ndisp - 1)/2 elements on each side.
"""
function vec2str(x :: AbstractVector{T}; ndisp :: Int=7) where T <: Union{AbstractFloat, Missing}
n = length(x)
if n ndisp
ndisp = n
nside = n
else
nside = max(1, div(ndisp - 1, 2))
end
s = "["
i = 1
while i nside
if x[i] !== missing
s *= @sprintf("%8.1e ", x[i])
else
s *= " ✗✗✗✗ "
end
i += 1
end
if i div(n, 2)
s *= "... "
end
i = max(i, n - nside + 1)
while i n
if x[i] !== missing
s *= @sprintf("%8.1e ", x[i])
else
s *= " ✗✗✗✗ "
end
i += 1
end
s *= "]"
return s
end
288 changes: 171 additions & 117 deletions test/test_aux.jl
Original file line number Diff line number Diff line change
@@ -1,119 +1,173 @@
@testset "aux" begin
# test Givens reflector corner cases
(c, s, ρ) = Krylov.sym_givens(0.0, 0.0)
@test (c == 1.0) && (s == 0.0) &&== 0.0)

a = 3.14
(c, s, ρ) = Krylov.sym_givens(a, 0.0)
@test (c == 1.0) && (s == 0.0) &&== a)
(c, s, ρ) = Krylov.sym_givens(-a, 0.0)
@test (c == -1.0) && (s == 0.0) &&== a)

b = 3.14
(c, s, ρ) = Krylov.sym_givens(0.0, b)
@test (c == 0.0) && (s == 1.0) &&== b)
(c, s, ρ) = Krylov.sym_givens(0.0, -b)
@test (c == 0.0) && (s == -1.0) &&== b)

(c, s, ρ) = Krylov.sym_givens(Complex(0.0), Complex(0.0))
@test (c == 1.0) && (s == Complex(0.0)) &&== Complex(0.0))

a = Complex(1.0, 1.0)
(c, s, ρ) = Krylov.sym_givens(a, Complex(0.0))
@test (c == 1.0) && (s == Complex(0.0)) &&== a)
(c, s, ρ) = Krylov.sym_givens(-a, Complex(0.0))
@test (c == 1.0) && (s == Complex(0.0)) &&== -a)

b = Complex(1.0, 1.0)
(c, s, ρ) = Krylov.sym_givens(Complex(0.0), b)
@test (c == 0.0) && (s == Complex(1.0)) &&== b)
(c, s, ρ) = Krylov.sym_givens(Complex(0.0), -b)
@test (c == 0.0) && (s == Complex(1.0)) &&== -b)

# test roots of a quadratic
roots = Krylov.roots_quadratic(0.0, 0.0, 0.0)
@test length(roots) == 1
@test roots[1] == 0.0

roots = Krylov.roots_quadratic(0.0, 0.0, 1.0)
@test length(roots) == 0

roots = Krylov.roots_quadratic(0.0, 3.14, -1.0)
@test length(roots) == 1
@test roots[1] == 1.0 / 3.14

roots = Krylov.roots_quadratic(1.0, 0.0, 1.0)
@test length(roots) == 0

roots = Krylov.roots_quadratic(1.0, 0.0, 0.0)
@test length(roots) == 2
@test roots[1] == 0.0
@test roots[2] == 0.0

roots = Krylov.roots_quadratic(1.0, 3.0, 2.0)
@test length(roots) == 2
@test roots[1] -2.0
@test roots[2] -1.0

roots = Krylov.roots_quadratic(1.0e+8, 1.0, 1.0)
@test length(roots) == 0

# ill-conditioned quadratic
roots = Krylov.roots_quadratic(-1.0e-8, 1.0e+5, 1.0, nitref=0)
@test length(roots) == 2
@test roots[1] == 1.0e+13
@test roots[2] == 0.0

# iterative refinement is crucial!
roots = Krylov.roots_quadratic(-1.0e-8, 1.0e+5, 1.0, nitref=1)
@test length(roots) == 2
@test roots[1] == 1.0e+13
@test roots[2] == -1.0e-05

# not ill-conditioned quadratic
roots = Krylov.roots_quadratic(-1.0e-7, 1.0, 1.0, nitref=0)
@test length(roots) == 2
@test isapprox(roots[1], 1.0e+7, rtol=1.0e-6)
@test isapprox(roots[2], -1.0, rtol=1.0e-6)

roots = Krylov.roots_quadratic(-1.0e-7, 1.0, 1.0, nitref=1)
@test length(roots) == 2
@test isapprox(roots[1], 1.0e+7, rtol=1.0e-6)
@test isapprox(roots[2], -1.0, rtol=1.0e-6)

# test trust-region boundary
x = ones(5)
d = ones(5); d[1:2:5] .= -1
@test_throws ErrorException Krylov.to_boundary(x, d, -1.0)
@test_throws ErrorException Krylov.to_boundary(x, d, 0.5)
@test_throws ErrorException Krylov.to_boundary(x, zeros(5), 1.0)
@test maximum(Krylov.to_boundary(x, d, 5.0)) 2.209975124224178
@test minimum(Krylov.to_boundary(x, d, 5.0)) -1.8099751242241782
@test maximum(Krylov.to_boundary(x, d, 5.0, flip=true)) 1.8099751242241782
@test minimum(Krylov.to_boundary(x, d, 5.0, flip=true)) -2.209975124224178

# test kzeros and kones
@test Krylov.kzeros(Vector{Float64}, 10) == zeros(10)
@test Krylov.kones(Vector{Float64}, 10) == ones(10)

# test ktypeof
a = rand(Float32, 10)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float32}
@test Krylov.ktypeof(b) == Vector{Float32}

a = rand(Float64, 10)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float64}
@test Krylov.ktypeof(b) == Vector{Float64}

a = sprand(Float32, 10, 0.5)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float32}
@test Krylov.ktypeof(b) == Vector{Float32}

a = sprand(Float64, 10, 0.5)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float64}
@test Krylov.ktypeof(b) == Vector{Float64}

@testset "sym_givens" begin
# test Givens reflector corner cases
(c, s, ρ) = Krylov.sym_givens(0.0, 0.0)
@test (c == 1.0) && (s == 0.0) &&== 0.0)

a = 3.14
(c, s, ρ) = Krylov.sym_givens(a, 0.0)
@test (c == 1.0) && (s == 0.0) &&== a)
(c, s, ρ) = Krylov.sym_givens(-a, 0.0)
@test (c == -1.0) && (s == 0.0) &&== a)

b = 3.14
(c, s, ρ) = Krylov.sym_givens(0.0, b)
@test (c == 0.0) && (s == 1.0) &&== b)
(c, s, ρ) = Krylov.sym_givens(0.0, -b)
@test (c == 0.0) && (s == -1.0) &&== b)

(c, s, ρ) = Krylov.sym_givens(Complex(0.0), Complex(0.0))
@test (c == 1.0) && (s == Complex(0.0)) &&== Complex(0.0))

a = Complex(1.0, 1.0)
(c, s, ρ) = Krylov.sym_givens(a, Complex(0.0))
@test (c == 1.0) && (s == Complex(0.0)) &&== a)
(c, s, ρ) = Krylov.sym_givens(-a, Complex(0.0))
@test (c == 1.0) && (s == Complex(0.0)) &&== -a)

b = Complex(1.0, 1.0)
(c, s, ρ) = Krylov.sym_givens(Complex(0.0), b)
@test (c == 0.0) && (s == Complex(1.0)) &&== b)
(c, s, ρ) = Krylov.sym_givens(Complex(0.0), -b)
@test (c == 0.0) && (s == Complex(1.0)) &&== -b)
end

@testset "roots_quadratic" begin
# test roots of a quadratic
roots = Krylov.roots_quadratic(0.0, 0.0, 0.0)
@test length(roots) == 1
@test roots[1] == 0.0

roots = Krylov.roots_quadratic(0.0, 0.0, 1.0)
@test length(roots) == 0

roots = Krylov.roots_quadratic(0.0, 3.14, -1.0)
@test length(roots) == 1
@test roots[1] == 1.0 / 3.14

roots = Krylov.roots_quadratic(1.0, 0.0, 1.0)
@test length(roots) == 0

roots = Krylov.roots_quadratic(1.0, 0.0, 0.0)
@test length(roots) == 2
@test roots[1] == 0.0
@test roots[2] == 0.0

roots = Krylov.roots_quadratic(1.0, 3.0, 2.0)
@test length(roots) == 2
@test roots[1] -2.0
@test roots[2] -1.0

roots = Krylov.roots_quadratic(1.0e+8, 1.0, 1.0)
@test length(roots) == 0

# ill-conditioned quadratic
roots = Krylov.roots_quadratic(-1.0e-8, 1.0e+5, 1.0, nitref=0)
@test length(roots) == 2
@test roots[1] == 1.0e+13
@test roots[2] == 0.0

# iterative refinement is crucial!
roots = Krylov.roots_quadratic(-1.0e-8, 1.0e+5, 1.0, nitref=1)
@test length(roots) == 2
@test roots[1] == 1.0e+13
@test roots[2] == -1.0e-05

# not ill-conditioned quadratic
roots = Krylov.roots_quadratic(-1.0e-7, 1.0, 1.0, nitref=0)
@test length(roots) == 2
@test isapprox(roots[1], 1.0e+7, rtol=1.0e-6)
@test isapprox(roots[2], -1.0, rtol=1.0e-6)

roots = Krylov.roots_quadratic(-1.0e-7, 1.0, 1.0, nitref=1)
@test length(roots) == 2
@test isapprox(roots[1], 1.0e+7, rtol=1.0e-6)
@test isapprox(roots[2], -1.0, rtol=1.0e-6)
end

@testset "to_boundary" begin
# test trust-region boundary
x = ones(5)
d = ones(5); d[1:2:5] .= -1
@test_throws ErrorException Krylov.to_boundary(x, d, -1.0)
@test_throws ErrorException Krylov.to_boundary(x, d, 0.5)
@test_throws ErrorException Krylov.to_boundary(x, zeros(5), 1.0)
@test maximum(Krylov.to_boundary(x, d, 5.0)) 2.209975124224178
@test minimum(Krylov.to_boundary(x, d, 5.0)) -1.8099751242241782
@test maximum(Krylov.to_boundary(x, d, 5.0, flip=true)) 1.8099751242241782
@test minimum(Krylov.to_boundary(x, d, 5.0, flip=true)) -2.209975124224178
end

@testset "kzeros" begin
# test kzeros
@test Krylov.kzeros(Vector{Float64}, 10) == zeros(Float64, 10)
@test Krylov.kzeros(Vector{ComplexF32}, 10) == zeros(ComplexF32, 10)
end

@testset "kones" begin
# test kones
@test Krylov.kones(Vector{Float64}, 10) == ones(Float64, 10)
@test Krylov.kones(Vector{ComplexF32}, 10) == ones(ComplexF32, 10)
end

@testset "ktypeof" begin
# test ktypeof
a = rand(Float32, 10)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float32}
@test Krylov.ktypeof(b) == Vector{Float32}

a = rand(Float64, 10)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float64}
@test Krylov.ktypeof(b) == Vector{Float64}

a = sprand(Float32, 10, 0.5)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float32}
@test Krylov.ktypeof(b) == Vector{Float32}

a = sprand(Float64, 10, 0.5)
b = view(a, 4:8)
@test Krylov.ktypeof(a) == Vector{Float64}
@test Krylov.ktypeof(b) == Vector{Float64}
end

@testset "macros" begin
# test macros
for FC (Float16, Float32, Float64, Complex{Float16}, Complex{Float32}, Complex{Float64})
n = 10
x = rand(FC, n)
y = rand(FC, n)
a = rand(FC)
b = rand(FC)
c = rand(FC)
s = rand(FC)

T = real(FC)
a2 = rand(T)
b2 = rand(T)

Krylov.@kdot(n, x, y)

Krylov.@kdotr(n, x, y)

Krylov.@knrm2(n, x)

Krylov.@kaxpy!(n, a, x, y)
Krylov.@kaxpy!(n, a2, x, y)

Krylov.@kaxpby!(n, a, x, b, y)
Krylov.@kaxpby!(n, a2, x, b, y)
Krylov.@kaxpby!(n, a, x, b2, y)
Krylov.@kaxpby!(n, a2, x, b2, y)

Krylov.@kcopy!(n, x, y)

Krylov.@kswap(x, y)

Krylov.@kref!(n, x, y, c, s)
end
end
end

0 comments on commit fc8677b

Please sign in to comment.