From d118cbd2b2cce7058c0bb0b57d7480c73bc2701e Mon Sep 17 00:00:00 2001 From: Pablo Winant Date: Thu, 28 Nov 2024 18:01:12 +0100 Subject: [PATCH] WIP: cleanup for DoloYAML --- Project.toml | 1 - src/adapt.jl | 14 ++++++------ src/dev_L2.jl | 31 ++++++++----------------- src/dolo_model.jl | 6 +++-- src/funs.jl | 5 ---- src/processes.jl | 58 ++++++++++++++++++++++++++++++++++++++++++++++- src/space.jl | 14 +++++++++++- 7 files changed, 90 insertions(+), 39 deletions(-) diff --git a/Project.toml b/Project.toml index 6abd21c..c532672 100644 --- a/Project.toml +++ b/Project.toml @@ -11,7 +11,6 @@ Crayons = "a8cc5b0e-0ffa-5ad4-8c14-923d3ee1735f" Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" Dolang = "e5c7262c-e9d2-5620-ad8e-1af14eb8a8e3" -DoloYAML = "5444aaa2-d204-4231-b785-bc52d56d058b" FiniteDiff = "6a86dc24-6348-571c-b903-95158fe2bd41" FiniteDifferences = "26cc04aa-876d-5657-8c51-4c34ba976000" Format = "1fa38f19-a742-5d3f-a2b9-30dd87b9d5f8" diff --git a/src/adapt.jl b/src/adapt.jl index 39873b6..d3f9f4d 100644 --- a/src/adapt.jl +++ b/src/adapt.jl @@ -14,13 +14,13 @@ function adapt_structure(to, v::GVector{G,V}) where G where V ) end -function adapt_structure(to, v::GArray{G,V}) where G where V - data = adapt(to, v.data) - GArray( - v.grid, - data - ) -end +# function adapt_structure(to, v::GArray{G,V}) where G where V +# data = adapt(to, v.data) +# GArray( +# v.grid, +# data +# ) +# end function adapt_structure(to, f::DFun{Dom, Gar, Itp, vars}) where Dom where Gar where Itp where vars itp = adapt(to, f.itp) values = adapt(to, f.values) diff --git a/src/dev_L2.jl b/src/dev_L2.jl index b61b11b..d1967ff 100644 --- a/src/dev_L2.jl +++ b/src/dev_L2.jl @@ -1,4 +1,8 @@ import Base: * +using LinearMaps +using BlockDiagonals +using LinearAlgebra: UniformScaling + struct LL{G,D_,F} grid::G @@ -10,29 +14,13 @@ struct I_L{_L} L::_L end - -using LinearMaps -using BlockDiagonals -using LinearAlgebra: UniformScaling - -(::UniformScaling, L::LL) = I_L(L) *(A::I_L, b) = b-A.L*b -convert(::Type{LinearMap}, IL::I_L) = let - I - convert(LinearMap,IL.L) - # L = IL.L - # # elt = eltype(L.D[1][1].F_x) - # elt = typeof(L.D[1][1].F_x) - # p,q = size(elt) - # N = length(L.grid) - # typ = SVector{q, Float64} - # fun = u->u-(ravel(L*GArray(L.grid, reinterpret(typ, u)))) - # LinearMap( - # fun, - # p*N, - # q*N - # ) -end +convert(::Type{LinearMap}, IL::I_L) = I - convert(LinearMap,IL.L) + + + # function *(A::I_L, b::AbstractVector) # x = GVector( @@ -49,9 +37,8 @@ end const LF = LL -import Base: ldiv! import Base: /,* -import LinearAlgebra: mul! +import LinearAlgebra: mul!, ldiv! # function mul!(L::LL, x::Number) diff --git a/src/dolo_model.jl b/src/dolo_model.jl index 4b6276f..89cfe2e 100644 --- a/src/dolo_model.jl +++ b/src/dolo_model.jl @@ -92,14 +92,16 @@ function discretize(model::AModel{<:VAR1}, d=Dict()) endo = get(d, :endo, Dict()) dvar = discretize(model.exogenous, exo) d = size(model.exogenous.Σ,1) + + # number of endogenous states n_s = length(Dolo.variables(model.states)) - d exo_grid = SGrid(dvar.Q) # TODO: simplify that call Tf = eltype(model) endo_space = CartesianSpace{n_s, Dolo.variables(model.states)[d+1:end],Tf}( - model.states.min[d+1:end], - model.states.max[d+1:end] + model.states.endo.min, + model.states.endo.max ) # return endo_space endo_grid = discretize(endo_space, endo) diff --git a/src/funs.jl b/src/funs.jl index f909a7b..e14af83 100644 --- a/src/funs.jl +++ b/src/funs.jl @@ -127,11 +127,6 @@ function (f::DFun{A,B,I,vars})(i::Int, j::Int) where A where B<:GArray{G,V} whe f.values[i,j] end -function (f::DFun{A,B,I,vars})(x::QP) where A where B<:GArray{G,V} where V where I where G<:PGrid{G1,G2} where G1<:SGrid where G2<:CGrid where vars - f(x.loc...) - # f(x.loc) -end - function (f::DFun{A,B,I,vars})(x::Tuple) where A where B<:GArray{G,V} where V where I where G<:PGrid{G1,G2} where G1<:SGrid where G2<:CGrid where vars f(x...) end diff --git a/src/processes.jl b/src/processes.jl index 02624ec..4a1aace 100644 --- a/src/processes.jl +++ b/src/processes.jl @@ -192,6 +192,10 @@ function convert_precision(T, var::Dolo.MvNormal{names,n,n2}) where n where n2 w end MvNormal(names, μ::SVector{n,T}, Σ::SMatrix{n,n, T,n2}) where n where n2 where T = MvNormal{names, n, n2,T}(μ,Σ) + +#special one-dimensional case +UNormal(names; μ=0.0, σ=1.0) = MvNormal(names, SVector(μ), SMatrix{1,1,typeof(μ),1}(σ^2)) + variables(mv::MvNormal{n}) where n = n MvNormal(names, Σ::SMatrix{n,n, T,n2}) where n where n2 where T = MvNormal{names,n,n2,T}(zero(SVector{n,T}),Σ) @@ -208,6 +212,11 @@ function MvNormal(names, Σ::Matrix) end +function MvNormal(names; Σ::Matrix=zeros(len(vars),len(vars))) + MvNormal(names, zeros(size(Σ, 1)), Σ) +end + + import Base: rand import Distributions @@ -248,7 +257,17 @@ struct VAR1{names,V,B} Σ::B end -VAR1(names, ρ, Σ) = VAR1{names, typeof(ρ), typeof(Σ)}(ρ, Σ) +VAR1(names; ρ=0.0, Σ=one(SMatrix{length(names),length(names),Float64,length(names)^2})) = VAR1(names,ρ,Σ) + +function VAR1(names, ρ, Σ) + Tf = typeof(ρ) + d = size(Σ,1) + + MT = SMatrix{d,d,Tf,d*d} + + VAR1{names, Tf, MT}(ρ,convert(MT,Σ)) +end + function rand(var::VAR1{N,V,B}, m::SVector{d,T}) where N where V where T where B<:SMatrix{1,1,Float64,1} where d @@ -344,6 +363,10 @@ ndims(mc::MarkovChain{names}) where names = length(names) # MarkovChain(names, P, Q) = MarkovChain{names, typeof(P), typeof(Q)}(P,Q) +function MarkovChain(names; transitions::Matrix, values::Matrix) + MarkovChain(names, transitions, values) +end + function MarkovChain(names, P::Matrix, Q::Matrix) Tf = eltype(P) d = size(P,1) @@ -364,3 +387,36 @@ MarkovProduct(mc::MarkovChain) = mc discretize(mc::MarkovChain, args...; kwargs...) = mc +# domains + + +struct ConstantProcess{names,C} + μ::C +end + +ConstantProcess(names; μ=zero(SVector{length(names),Float64})) = ConstantProcess{names,typeof(μ)}(SVector(μ...)) + +variables(cp::ConstantProcess{vars}) where vars = vars + +function get_domain(cp::ConstantProcess) + d = length(cp.μ) + T = eltype(cp.μ) + min = fill(-T(Inf), d) ### not absolutely clear this is the right default, certainly not! + max = fill(T(Inf), d) + names = variables(cp) + return Dolo.CSpace(; (k=>(a,b) for (k,a,b) in zip(names,min,max))...) +end + +function get_domain(iid::P) where P<:Union{VAR1,MvNormal} + d = size(iid.Σ,1) + T = eltype(iid.Σ) + min = fill(-T(Inf), d) ### not absolutely clear this is the right default, certainly not! + max = fill(T(Inf), d) + names = variables(iid) + return Dolo.CSpace(; (k=>(a,b) for (k,a,b) in zip(names,min,max))...) +end + + +function get_domain(mc::MarkovChain) + GridSpace(variables(mc),mc.Q) +end diff --git a/src/space.jl b/src/space.jl index 34666da..3745a1a 100644 --- a/src/space.jl +++ b/src/space.jl @@ -8,6 +8,8 @@ end const CSpace = CartesianSpace +ndims(::CartesianSpace{d}) where d = d + # CartesianSpace(a::NTuple{d,Tf}, b::NTuple{d, Tf}) where d where Tf = CartesianSpace{length(a), (:x,), Tf}(a,b) function CartesianSpace(kwargs::Pair{Symbol, Tuple{Tf, Tf}}...) where Tf @@ -41,7 +43,7 @@ Base.in(e::SVector, cs) = all( ( (e[i]<=cs.max[i])&(e[i]>=cs.min[i]) for i=1:len # dims(dom::CartesianSpace{d,dims}) where d where dims = dims -ndims(dom::CartesianSpace{d, dims}) where d where dims = d +# ndims(dom::CartesianSpace{d, dims}) where d where dims = d variables(dom::CartesianSpace{d,t}) where d where t = t dims(dom::CartesianSpace) = variables(dom) @@ -75,7 +77,17 @@ end ProductSpace(A,B) = ProductSpace((A,B)) +import Base +function Base.getproperty(PS::ProductSpace,sym::Symbol) + if sym==:exo + return PS.spaces[1] + elseif sym==:endo + return PS.spaces[2] + else + getfield(PS, sym) + end +end function draw(p::ProductSpace) a = rand(p.spaces[1])