Skip to content

Commit

Permalink
works but not optimal
Browse files Browse the repository at this point in the history
  • Loading branch information
@thorek1
thorek1 committed Nov 10, 2023
1 parent fda69ed commit 3112064
Show file tree
Hide file tree
Showing 3 changed files with 93 additions and 46 deletions.
65 changes: 53 additions & 12 deletions src/MacroModelling.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -123,6 +123,9 @@ Base.show(io::IO, 𝓂::β„³) = println(io,
# "\nNon-stochastic-steady-state found: ",!𝓂.solution.outdated_NSSS
)

check_for_dynamic_variables(ex::Int) = false
check_for_dynamic_variables(ex::Float64) = false
check_for_dynamic_variables(ex::Symbol) = occursin(r"β‚β‚β‚Ž|β‚β‚€β‚Ž|β‚β‚‹β‚β‚Ž",string(ex))

function check_for_dynamic_variables(ex::Expr)
dynamic_indicator = Bool[]
Expand Down Expand Up @@ -247,7 +250,12 @@ function set_up_obc_violation_function!(𝓂)
$(steady_state...)
$(steady_state_obc...)

return_value = $(𝓂.obc_violation_equations...)
constraint_values = Vector{JuMP.AffExpr}[]
shock_sign_indicators = Bool[]

$(𝓂.obc_violation_equations...)

return vcat(constraint_values...), shock_sign_indicators
end)

𝓂.obc_violation_function = @RuntimeGeneratedFunction(calc_obc_violation)
Expand Down Expand Up @@ -289,20 +297,51 @@ function write_obc_violation_equations(𝓂)
begin
plchldr = Symbol(replace(string(x.args[2]), "β‚β‚€β‚Ž" => ""))

ineq_plchldr_1 = Symbol(replace(string(x.args[3].args[2]), "β‚β‚€β‚Ž" => ""))
ineq_plchldr_1 = x.args[3].args[2] isa Symbol ? Symbol(replace(string(x.args[3].args[2]), "β‚β‚€β‚Ž" => "")) : x.args[3].args[2]

:(if isapprox($plchldr, $ineq_plchldr_1, atol = 1e-12)
if contains($(string(plchldr)), "⁺")
return $(x.args[3].args[3]), false
else
return $(x.args[3].args[2]), true
arg1 = x.args[3].args[2]
arg2 = x.args[3].args[3]

dyn_1 = check_for_dynamic_variables(x.args[3].args[2])
dyn_2 = check_for_dynamic_variables(x.args[3].args[3])

cond1 = Expr[]
cond2 = Expr[]

maximisation = contains(string(plchldr), "⁺")

if dyn_1
if maximisation
push!(cond1, :(push!(constraint_values, $(x.args[3].args[2]))))
push!(cond2, :(push!(constraint_values, $(x.args[3].args[2]))))
else
push!(cond1, :(push!(constraint_values, -$(x.args[3].args[2]))))
push!(cond2, :(push!(constraint_values, -$(x.args[3].args[2])))) # RBC
end
else
if contains($(string(plchldr)), "⁻")
return -$(x.args[3].args[2]), true
else
return -$(x.args[3].args[3]), false
end

if dyn_2
if maximisation
push!(cond1, :(push!(constraint_values, $(x.args[3].args[3]))))
push!(cond2, :(push!(constraint_values, $(x.args[3].args[3])))) # testmax
else
push!(cond1, :(push!(constraint_values, -$(x.args[3].args[3]))))
push!(cond2, :(push!(constraint_values, -$(x.args[3].args[3])))) # RBC
end
end

if maximisation
push!(cond1, :(push!(shock_sign_indicators, true)))
push!(cond2, :(push!(shock_sign_indicators, true)))
else
push!(cond1, :(push!(shock_sign_indicators, false)))
push!(cond2, :(push!(shock_sign_indicators, false)))
end

:(if isapprox($plchldr, $ineq_plchldr_1, atol = 1e-12)
$(Expr(:block, cond1...))
else
$(Expr(:block, cond2...))
end)
end :
x :
Expand Down Expand Up @@ -3508,6 +3547,8 @@ write_parameters_input!(𝓂::β„³, parameters::Pair{String,Float64}; verbose::Bo


write_parameters_input!(𝓂::β„³, parameters::Tuple{Pair{Symbol,Float64},Vararg{Pair{Symbol,Float64}}}; verbose::Bool = true) = write_parameters_input!(𝓂::β„³, Dict(parameters), verbose = verbose)
write_parameters_input!(𝓂::β„³, parameters::Tuple{Pair{Symbol,Float64},Vararg{Pair{String,Int}}}; verbose::Bool = true) = write_parameters_input!(𝓂::β„³, Dict(parameters), verbose = verbose)
write_parameters_input!(𝓂::β„³, parameters::Tuple{Pair{Symbol,Int},Vararg{Pair{String,Float64}}}; verbose::Bool = true) = write_parameters_input!(𝓂::β„³, Dict(parameters), verbose = verbose)
write_parameters_input!(𝓂::β„³, parameters::Tuple{Pair{String,Float64},Vararg{Pair{String,Float64}}}; verbose::Bool = true) = write_parameters_input!(𝓂::β„³, Dict([i[1] |> Meta.parse |> replace_indices => i[2] for i in parameters])
, verbose = verbose)
write_parameters_input!(𝓂::β„³, parameters::Vector{Pair{Symbol, Float64}}; verbose::Bool = true) = write_parameters_input!(𝓂::β„³, Dict(parameters), verbose = verbose)
Expand Down
37 changes: 20 additions & 17 deletions src/get_functions.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -912,7 +912,8 @@ function get_irf(𝓂::β„³;
# Find shocks fulfilling constraint
# model = JuMP.Model(MadNLP.Optimizer)
model = JuMP.Model(NLopt.Optimizer)
JuMP.set_attribute(model, "algorithm", :LD_MMA)
JuMP.set_attribute(model, "algorithm", :LD_SLSQP)
# JuMP.set_attribute(model, "algorithm", :LD_MMA)

JuMP.set_silent(model)

Expand All @@ -922,22 +923,24 @@ function get_irf(𝓂::β„³;
JuMP.@variable(model, x[1:num_shocks*periods_per_shock])

# Now loop through obc_shock_bounds to set the bounds on these variables.
for (idx, v) in enumerate(𝓂.var[obc_inequalities_idx])
idxs = (idx - 1) * periods_per_shock + 1:idx * periods_per_shock
if contains(string(v), "α΅’α΅‡αΆœβΊ")
if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2]
JuMP.set_upper_bound.(x[idxs], 0)
else
JuMP.set_lower_bound.(x[idxs], 0)
end
else
if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2]
JuMP.set_lower_bound.(x[idxs], 0)
else
JuMP.set_upper_bound.(x[idxs], 0)
end
end
end
# for (idx, v) in enumerate(𝓂.var[obc_inequalities_idx])
# idxs = (idx - 1) * periods_per_shock + 1:idx * periods_per_shock
# # if contains(string(v), "α΅’α΅‡αΆœβΊ")
# if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2][idx]
# # JuMP.set_upper_bound.(x[idxs], 0)
# JuMP.set_lower_bound.(x[idxs], 0)
# else
# JuMP.set_upper_bound.(x[idxs], 0)
# # JuMP.set_lower_bound.(x[idxs], 0)
# end
# # else
# # if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2][idx]
# # JuMP.set_lower_bound.(x[idxs], 0)
# # else
# # JuMP.set_upper_bound.(x[idxs], 0)
# # end
# # end
# end

JuMP.@objective(model, Min, x' * β„’.I * x)

Expand Down
37 changes: 20 additions & 17 deletions src/plotting.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -473,7 +473,8 @@ function plot_irf(𝓂::β„³;
# Find shocks fulfilling constraint
# model = JuMP.Model(MadNLP.Optimizer)
model = JuMP.Model(NLopt.Optimizer)
JuMP.set_attribute(model, "algorithm", :LD_MMA)
JuMP.set_attribute(model, "algorithm", :LD_SLSQP)
# JuMP.set_attribute(model, "algorithm", :LD_MMA)

JuMP.set_silent(model)

Expand All @@ -483,22 +484,24 @@ function plot_irf(𝓂::β„³;
JuMP.@variable(model, x[1:num_shocks*periods_per_shock])

# Now loop through obc_shock_bounds to set the bounds on these variables.
for (idx, v) in enumerate(𝓂.var[obc_inequalities_idx])
idxs = (idx - 1) * periods_per_shock + 1:idx * periods_per_shock
if contains(string(v), "α΅’α΅‡αΆœβΊ")
if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2]
JuMP.set_upper_bound.(x[idxs], 0)
else
JuMP.set_lower_bound.(x[idxs], 0)
end
else
if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2]
JuMP.set_lower_bound.(x[idxs], 0)
else
JuMP.set_upper_bound.(x[idxs], 0)
end
end
end
# for (idx, v) in enumerate(𝓂.var[obc_inequalities_idx])
# idxs = (idx - 1) * periods_per_shock + 1:idx * periods_per_shock
# # if contains(string(v), "α΅’α΅‡αΆœβΊ")
# if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2][idx]
# # JuMP.set_upper_bound.(x[idxs], 0)
# JuMP.set_lower_bound.(x[idxs], 0)
# else
# JuMP.set_upper_bound.(x[idxs], 0)
# # JuMP.set_lower_bound.(x[idxs], 0)
# end
# # else
# # if 𝓂.obc_violation_function(x, past_states, past_shocks, state_update, reference_steady_state, 𝓂, unconditional_forecast_horizon, JuMP.AffExpr.(present_shocks))[2][idx]
# # JuMP.set_lower_bound.(x[idxs], 0)
# # else
# # JuMP.set_upper_bound.(x[idxs], 0)
# # end
# # end
# end

JuMP.@objective(model, Min, x' * β„’.I * x)

Expand Down

0 comments on commit 3112064

Please sign in to comment.