QuantEcon / BasisMatrices.jl / 147

Committed 6 Dec 2017 - 21:20 coverage: 92.635%. First build

Build # 147

Build Type

Pull #46

travis-ci

Committed by

web-flow

Commit Message

ENH: full/proper support for multi-dimensional BasisParams

closes #45

Pull Request Pull Request #46: ENH: full/proper support for multi-dimensional BasisParams

Run Details

124 of 132 new or added lines in 7 files covered. (93.94%)

1283 of 1385 relevant lines covered (92.64%)

506206.2 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

87.18

/src/smolyak.jl

#=

Chase Coleman (ccoleman@stern.nyu.edu)
Spencer Lyon (spencer.lyon@stern.nyu.edu)

References
----------

Judd, Kenneth L, Lilia Maliar, Serguei Maliar, and Rafael Valero. 2013.
    "Smolyak Method for Solving Dynamic Economic Models: Lagrange
    Interpolation, Anisotropic Grid and Adaptive Domain".

Krueger, Dirk, and Felix Kubler. 2004. "Computing Equilibrium in OLG
    Models with Stochastic Production." Journal of Economic Dynamics and
    Control 28 (7) (April): 1411-1436.

=#

# include helper file
include("smol_util.jl")

struct Smolyak <: BasisFamily end

struct SmolyakParams{T,Tmu<:IntSorV,d} <: BasisParams
    mu::Tmu
    a::Vector{T}
    b::Vector{T}

    # internal fields
    inds::Vector{Vector{Int}}  # Smolyak indices
    pinds::Matrix{Int64}  # Polynomial indices

    function SmolyakParams{T,Tmu,d}(
            mu::Tmu, a::Vector{T}, b::Vector{T}
        ) where {T,Tmu,d}
        d < 2 && error("You passed d = $d. d must be greater than 1")
        if length(mu) > 1
            # working on building an anisotropic grid
            if length(mu) != d
                error("mu must have d elements. It has $(length(mu))")
            end

            if any(mu .< 1)
                error("All elements of mu must be greater than 1 (recieved $mu)")
            end
        else
            mu < 1 && error("You passed mu = $mu. mu must be greater than 1")
        end

        inds = smol_inds(d, mu)
        pinds = poly_inds(d, mu, inds)

        new{T,Tmu,d}(mu, a, b, inds, pinds)
    end
end

function SmolyakParams(
        ::Union{Val{d},Type{Val{d}}}, mu::Tmu, a::Vector{T}=fill(-1.0, d),
        b::Vector{T}=fill(1.0, d)
    ) where {T,Tmu,d}
    SmolyakParams{T,Tmu,d}(mu, a, b)
end

# add methods to helper routines from other files
smol_inds(sp::SmolyakParams) = sp.inds

function poly_inds(
        sp::SmolyakParams{T,Tmu,d}, inds::Vector{Vector{Int}}=sp.inds
    ) where {T,Tmu,d}
    poly_inds(d, sp.mu, inds)
end

function build_grid(
        sp::SmolyakParams{T,Tmu,d}, inds::Vector{Vector{Int}}=sp.inds
    ) where {T,Tmu,d}
    build_grid(d, sp.mu, inds)
end

function build_B!(
        out::AbstractMatrix{T}, sp::SmolyakParams{Tp,Tmu,d},
        pts::Matrix{Float64}, b_inds::Matrix{Int64}=sp.pinds
    ) where {T,Tp,Tmu,d}
    build_B!(out, d, sp.mu, pts, b_inds)
end

function build_B(sp::SmolyakParams, pts::Matrix{Float64},
                 b_inds::Matrix{Int64}=sp.pinds)
    build_B!(Array{Float64}(size(pts, 1), size(b_inds, 1)), sp, pts, b_inds)
end

function dom2cube!(out::AbstractMatrix, pts::AbstractMatrix,
                   sp::SmolyakParams)
    dom2cube!(out, pts, sp.a, sp.b)
end
dom2cube(pts::AbstractMatrix, sp::SmolyakParams) = dom2cube(pts, sp.a, sp.b)

function cube2dom!(out::AbstractMatrix, pts::AbstractMatrix,
                   sp::SmolyakParams)
    cube2dom!(out, pts, sp.a, sp.b)
end
cube2dom(pts::AbstractMatrix, sp::SmolyakParams) = cube2dom(pts, sp.a, sp.b)

## BasisParams interface
# define these methods on the type, the instance version is defined over
# BasisParams
family(::Type{T}) where {T<:SmolyakParams} = Smolyak
family_name(::Type{T}) where {T<:SmolyakParams} = "Smolyak"
Base.eltype(::Type{SmolyakParams{T,Tmu,d}}) where {T,Tmu,d} = T

# methods that only make sense for instances
Base.min(p::SmolyakParams) = p.a
Base.max(p::SmolyakParams) = p.b
Base.ndims(::Union{SmolyakParams{T,Tmu,d},Type{SmolyakParams{T,Tmu,d}}}) where {T,Tmu,d} = d

function Base.show(io::IO, p::SmolyakParams{T,Tmu,d}) where {T,Tmu,d}
    m = string("Smolyak interpolation parameters in $(d) dimensions",
               " from $(p.a) �� $(p.b)")
    print(io, m)
end

# TODO: fix this
function Base.length(sp::SmolyakParams{T,Tmu,d})::Int where {T,Tmu,d}
    mu = sp.mu
    mu == 1 ? 2d - 1 :
    mu == 2 ? Int(1 + 4d + 4d*(d-1)/2 ) :
    mu == 3 ? Int(1 + 8d + 12d*(d-1)/2 + 8d*(d-1)*(d-2)/6) :
    error("We only know the number of grid points for mu ��� [1, 2, 3]")
end

function nodes(sp::SmolyakParams)
    dom_grid = build_grid(sp)
    cube2dom!(dom_grid, dom_grid, sp)
end

function evalbase(sp::SmolyakParams, x::AbstractVector, order=0)
    evalbase(sp, reshape(x, 1, :), order)
end

function evalbase(sp::SmolyakParams, x::AbstractArray, order=0)
    for _o in order
        if _o != 0
        m = string(
            "Derivatives for Smolyak are not implemented yet...",
            "\nPlease open an issue on GitHub if you would like use this feature"
        )
        error(m)
        end
    end
    cube_pts = dom2cube(x, sp)
    build_B(sp, cube_pts, sp.pinds)
end

1	#=
2
3	Chase Coleman (ccoleman@stern.nyu.edu)
4	Spencer Lyon (spencer.lyon@stern.nyu.edu)
5
6	References
7	----------
8
9	Judd, Kenneth L, Lilia Maliar, Serguei Maliar, and Rafael Valero. 2013.
10	"Smolyak Method for Solving Dynamic Economic Models: Lagrange
11	Interpolation, Anisotropic Grid and Adaptive Domain".
12
13	Krueger, Dirk, and Felix Kubler. 2004. "Computing Equilibrium in OLG
14	Models with Stochastic Production." Journal of Economic Dynamics and
15	Control 28 (7) (April): 1411-1436.
16
17	=#
18
19	# include helper file
20	include("smol_util.jl")
21
22	struct Smolyak <: BasisFamily end
23
24	struct SmolyakParams{T,Tmu<:IntSorV,d} <: BasisParams
25	mu::Tmu
26	a::Vector{T}
27	b::Vector{T}
28
29	# internal fields
30	inds::Vector{Vector{Int}} # Smolyak indices
31	pinds::Matrix{Int64} # Polynomial indices
32
33	function SmolyakParams{T,Tmu,d}(
34	mu::Tmu, a::Vector{T}, b::Vector{T}
35	) where {T,Tmu,d}
36	d < 2 && error("You passed d = $d. d must be greater than 1")	10×
37	if length(mu) > 1	9×
38	# working on building an anisotropic grid
39	if length(mu) != d	4×
40	error("mu must have d elements. It has $(length(mu))")	1×
41	end
42
43	if any(mu .< 1)	3×
44	error("All elements of mu must be greater than 1 (recieved $mu)")	1×
45	end
46	else
47	mu < 1 && error("You passed mu = $mu. mu must be greater than 1")	5×
48	end
49
50	inds = smol_inds(d, mu)	6×
51	pinds = poly_inds(d, mu, inds)	6×
52
53	new{T,Tmu,d}(mu, a, b, inds, pinds)	6×
54	end
55	end
56
57	function SmolyakParams(
58	::Union{Val{d},Type{Val{d}}}, mu::Tmu, a::Vector{T}=fill(-1.0, d),
59	b::Vector{T}=fill(1.0, d)
60	) where {T,Tmu,d}
61	SmolyakParams{T,Tmu,d}(mu, a, b)	16×
62	end
63
64	# add methods to helper routines from other files
65	smol_inds(sp::SmolyakParams) = sp.inds	5×
66
67	function poly_inds(
68	sp::SmolyakParams{T,Tmu,d}, inds::Vector{Vector{Int}}=sp.inds
69	) where {T,Tmu,d}
70	poly_inds(d, sp.mu, inds)	10×
71	end
72
73	function build_grid(
74	sp::SmolyakParams{T,Tmu,d}, inds::Vector{Vector{Int}}=sp.inds
75	) where {T,Tmu,d}
76	build_grid(d, sp.mu, inds)	30×
77	end
78
79	function build_B!(
80	out::AbstractMatrix{T}, sp::SmolyakParams{Tp,Tmu,d},
81	pts::Matrix{Float64}, b_inds::Matrix{Int64}=sp.pinds
82	) where {T,Tp,Tmu,d}
83	build_B!(out, d, sp.mu, pts, b_inds)	29×
84	end
85
86	function build_B(sp::SmolyakParams, pts::Matrix{Float64},
87	b_inds::Matrix{Int64}=sp.pinds)
88	build_B!(Array{Float64}(size(pts, 1), size(b_inds, 1)), sp, pts, b_inds)	24×
89	end
90
91	function dom2cube!(out::AbstractMatrix, pts::AbstractMatrix,
92	sp::SmolyakParams)
93	dom2cube!(out, pts, sp.a, sp.b)	5×
94	end
95	dom2cube(pts::AbstractMatrix, sp::SmolyakParams) = dom2cube(pts, sp.a, sp.b)	24×
96
97	function cube2dom!(out::AbstractMatrix, pts::AbstractMatrix,
98	sp::SmolyakParams)
99	cube2dom!(out, pts, sp.a, sp.b)	20×
100	end
101	cube2dom(pts::AbstractMatrix, sp::SmolyakParams) = cube2dom(pts, sp.a, sp.b)	10×
102
103	## BasisParams interface
104	# define these methods on the type, the instance version is defined over
105	# BasisParams
106	family(::Type{T}) where {T<:SmolyakParams} = Smolyak	5×
107	family_name(::Type{T}) where {T<:SmolyakParams} = "Smolyak"	5×
108	Base.eltype(::Type{SmolyakParams{T,Tmu,d}}) where {T,Tmu,d} = T	7×
109
110	# methods that only make sense for instances
111	Base.min(p::SmolyakParams) = p.a	10×
112	Base.max(p::SmolyakParams) = p.b	10×
113	Base.ndims(::Union{SmolyakParams{T,Tmu,d},Type{SmolyakParams{T,Tmu,d}}}) where {T,Tmu,d} = d	140×
114
115	function Base.show(io::IO, p::SmolyakParams{T,Tmu,d}) where {T,Tmu,d}
116	m = string("Smolyak interpolation parameters in $(d) dimensions",	1×
117	" from $(p.a) �� $(p.b)")
118	print(io, m)	1×
119	end
120
121	# TODO: fix this
122	function Base.length(sp::SmolyakParams{T,Tmu,d})::Int where {T,Tmu,d}
NEW 123	mu = sp.mu	!
124	mu == 1 ? 2d - 1 :	!
125	mu == 2 ? Int(1 + 4d + 4d*(d-1)/2 ) :
126	mu == 3 ? Int(1 + 8d + 12d(d-1)/2 + 8d(d-1)*(d-2)/6) :
127	error("We only know the number of grid points for mu �� [1, 2, 3]")
128	end
129
130	function nodes(sp::SmolyakParams)
131	dom_grid = build_grid(sp)	15×
132	cube2dom!(dom_grid, dom_grid, sp)	15×
133	end
134
135	function evalbase(sp::SmolyakParams, x::AbstractVector, order=0)
136	evalbase(sp, reshape(x, 1, :), order)	!
137	end
138
139	function evalbase(sp::SmolyakParams, x::AbstractArray, order=0)
140	for _o in order	24×
141	if _o != 0	16×
NEW 142	m = string(	!
143	"Derivatives for Smolyak are not implemented yet...",
144	"\nPlease open an issue on GitHub if you would like use this feature"
145	)
NEW 146	error(m)	!
147	end
148	end
149	cube_pts = dom2cube(x, sp)	14×
150	build_B(sp, cube_pts, sp.pinds)	14×
151	end

QuantEcon / BasisMatrices.jl / 147

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous