#
#--------------------------------------#
export restrict, restrictUpper, restrictLower
#--------------------------------------#
"""
 R * u = u[2:end-1]
 R'* u = [0,u,0]
"""
function restrict(n::Integer)
    Id = sparse(I,n,n)
    R  = Id[2:n-1,:]
    return R
end
"""
 R * u = u[2:end]
 R'* u = [0,u]
"""
function restrictUpper(n::Integer)
    Id = sparse(I,n,n)
    R  = Id[2:n,:]
    return R
end
"""
 R * u = u[1:end-1]
 R'* u = [u,0]
"""
function restrictLower(n::Integer)
    Id = sparse(I,n,n)
    R  = Id[1:n-1,:]
    return R
end
#--------------------------------------#
export Tdg
#--------------------------------------#
"""
 Adjointable tridiagonal functor
"""
struct Tdg
    dl::AbstractVector # lower diag
    dd::AbstractVector # main  diag
    du::AbstractVector # upper diag
end

Flux.@functor Tdg

function (T::Tdg)(u::AbstractArray)
    n  = length(T.dd)

    Ru,Rl = Zygote.ignore() do
        Ru = restrictUpper(n)
        Rl = restrictLower(n)
        return Ru,Rl
    end

    v  = u .* T.dd              # main  diag
    v += Rl' * (T.du .* (Ru*u)) # upper diag
    v += Ru' * (T.dl .* (Rl*u)) # lower diag

    return v
end

Base.:*(T::Tdg,u::AbstractArray) = T(u)
Base.:adjoint(T::Tdg) = Tdg(T.du,T.dd,T.dl)
#--------------------------------------#
export lsq
#--------------------------------------#
function lsq(u,A)  # target data, linear function

    x = A \ u

#   AA = Matrix{Float32}(A)
#   uu = vec(u')
#   x,stat = Krylov.cgls(AA,uu)

    return x
end
#--------------------------------------#
#
