diff --git a/misc/dev_float32.jl b/misc/dev_float32.jl index e7f5456..f074c9d 100644 --- a/misc/dev_float32.jl +++ b/misc/dev_float32.jl @@ -1,8 +1,9 @@ using Dolo root_dir = pkgdir(Dolo) -model = include("$(root_dir)/examples/ymodels/rbc_mc.jl") +model32 = include("$(root_dir)/misc/rbc_float32.jl") -dm = Dolo.discretize(model, Dict(:endo=>[10000000]) ) + +dm32 = Dolo.discretize(model, Dict(:endo=>[10000000]) ) Dolo.convert \ No newline at end of file diff --git a/misc/rbc_float32.jl b/misc/rbc_float32.jl index 3958bab..d6beddf 100644 --- a/misc/rbc_float32.jl +++ b/misc/rbc_float32.jl @@ -56,22 +56,47 @@ model = let end -dmodel = Dolo.discretize(model) +# now these are orphans +model32 = Dolo.convert_precision(Float32,model) -model32 = Dolo.convert_precision(Float32,model) -dmodel32 = Dolo.discretize(model32) +function Dolo.transition(model::typeof(model32), s::NamedTuple, x::NamedTuple, M::NamedTuple) + + (;δ, ρ) = model.calibration + + # Z = e.Z + K = s.k * (1-δ) + x.i -wksp = Dolo.time_iteration_workspace(dmodel32) + (;k=K,) ## This is only the endogenous state -itps = wksp.φ.itp +end -(;x0,r0, φ) = wksp -Dolo.F(dmodel32, x0, φ) +function intermediate(model::typeof(model32),s::NamedTuple, x::NamedTuple) + + p = model.calibration + y = exp(s.z)*(s.k^p.α)*(x.n^(1-p.α)) + w = (1-p.α)*y/x.n + rk = p.α*y/s.k + c = y - x.i + return ( (; y, c, rk, w)) + +end +function arbitrage(model::typeof(model32), s::NamedTuple, x::NamedTuple, S::NamedTuple, X::NamedTuple) + + p = model.calibration + + y = intermediate(model, s, x) + Y = intermediate(model, S, X) + res_1 = p.χ*(x.n^p.η)*(y.c^p.σ) - y.w + res_2 = (p.β*(y.c/Y.c)^p.σ)*(1 - p.δ + Y.rk) - 1 + + return ( (;res_1, res_2) ) + +end -Dolo.time_iteration(model) \ No newline at end of file +model32 \ No newline at end of file diff --git a/src/algos/time_iteration.jl b/src/algos/time_iteration.jl index 3ecc715..c8ec90c 100644 --- a/src/algos/time_iteration.jl +++ b/src/algos/time_iteration.jl @@ -195,6 +195,7 @@ function time_iteration(model::DYModel, lam = 0.5 local η_0 = NaN + local η convergence = false iterations = T @@ -246,6 +247,7 @@ function time_iteration(model::DYModel, ε_n = norm(r0) if ε_nmaximum(abs, u), v.data) @@ -45,7 +45,10 @@ end eltype(g::GArray{G,T}) where G where T = eltype(T) # warning: these functions don't copy any data -ravel(g::GArray) = reinterpret(Float64, g.data) +function ravel(g::GArray) + Tf = eltype(g.data[1]) + reinterpret(Tf, g.data) +end unravel(g::GArray, x) = GArray( g.grid, reinterpret(eltype(g), x) @@ -113,22 +116,21 @@ import Base: *, \, +, -, / *(x::Number, A::GArray{G,T}) where G where T = GArray(A.grid, x .* A.data) -*(A::GArray{G,Vector{T}}, x::SVector{q, Float64}) where G where T <:SMatrix{p, q, Float64, n} where p where q where n = GArray(A.grid, [M*x for M in A.data]) -# *(A::GArray{G,Vector{T}}, x::SLArray{Tuple{q}, Float64, 1, q, names}) where G where T <:SMatrix{p, q, Float64, n} where p where q where n where names = A*SVector(x...) +*(A::GArray{G,Vector{T}}, x::SVector{q, Tf}) where G where T <:SMatrix{p, q, Tf, n} where p where q where n where Tf = GArray(A.grid, [M*x for M in A.data]) -*(A::GArray{G,T}, B::AbstractArray{Float64}) where G where T <:SMatrix{p, q, Float64, n} where p where q where n = +*(A::GArray{G,T}, B::AbstractArray{Tf}) where G where T <:SMatrix{p, q, Tf, n} where p where q where n where Tf = ravel( GArray( A.grid, - A.data .* reinterpret(SVector{q, Float64}, B) + A.data .* reinterpret(SVector{q, Tf}, B) ) ) import Base: convert -function Base.convert(::Type{Matrix}, A::GArray{G,Vector{T}}) where G where T <:SMatrix{p, q, Float64, k} where p where q where k +function Base.convert(::Type{Matrix}, A::GArray{G,Vector{T}}) where G where T <:SMatrix{p, q, Tf, k} where p where q where k where Tf N = length(A.data) n0 = N*p n1 = N*q diff --git a/src/grids.jl b/src/grids.jl index 3ccec25..07062c7 100644 --- a/src/grids.jl +++ b/src/grids.jl @@ -3,7 +3,7 @@ abstract type ASGrid{d} <: AGrid{d} end import Base: eltype, iterate, size -eltype(cg::AGrid{d}) where d = SVector{d, Float64} +# eltype(cg::AGrid{d}) where d = SVector{d, Float64} ndims(cg::AGrid{d}) where d = d struct CGrid{d,Tf} <: AGrid{d} @@ -90,7 +90,8 @@ from_linear(g::PGrid{G1, G2, d}, n) where G1 where G2 where d = let x=divrem(n-1 getindex(g::PGrid{G1, G2, d}, n::Int) where G1 where G2 where d = getindex(g, from_linear(g, n)...) function getindex(g::PGrid{G1, G2, d}, i::Int64, j::Int64) where G1<:SGrid{d1} where G2<:CGrid{d2} where d where d1 where d2 - SVector{d,Float64}(g.grids[1][i]..., g.grids[2][j]...) + Tf = eltype(g) + SVector{d,Tf}(g.grids[1][i]..., g.grids[2][j]...) end @@ -177,17 +178,21 @@ end function Base.iterate(g::PGrid{G1, G2, d}) where G1 where G2 where d + T = eltype(g) x = g.grids[1][1] y = g.grids[2][1] - return (SVector{d, Float64}(x...,y...),(y,1,1)) + return (SVector{d, T}(x...,y...),(y,1,1)) end function Base.iterate(g::PGrid{G1,G2,d},state) where G1 where G2 where d + + T = eltype(g) + y,i,j=state if i