# Updating Values via Iteration

I’m having a few issues with using a while loop to update values and return. The code seems to work fine when wrapped in functions but not otherwise. As I’m still learning the language I want to understand where I’m going wrong. I’m fairly sure it’s to do with scoping rules that I don’t fully understand.

The issue is that whenever I run the code it seems to converge at (update: It does indeed converge at iteration 2 with error 0.0 while the function based code converges at iteration 5 with error 0.0) iteration 1 with error 0.0 even though the pre and updated values are clearly different. Somewhere in my code I must be setting the pre values equal to the updated without intending too. Please see the code below. The variables of interest are policy_new and policy. I apologise if the code is rather poor; as I mentioned I’m still learning. Any help and advice is greatly appreciated. Feel free to comment on all aspects.

``````using Plots
using SparseArrays
using LinearAlgebra

#                     Model Parameters, Grid, Utility

β = 0.95;
α = 0.65;
δ = 1;
cover = 0.1;
N = 5;
k_ss = ((1/β-1+δ)/α)^(1/(α-1));
k_grid = range((1-cover)*k_ss, (1+cover)*k_ss, length=N);
max_iter = 500;
diff = Inf;
tol = 1e-6;
u(c) =  c <= 0 ? -Inf : log(c);

function find_maximiser(A)
v, B = findmax(A, dims=2)
B = getindex.(B,2)
return v, B
end

#                    Start Iterations
policy = k_grid
v = u.(k_grid.^α)./(1-β)

for iter ∈ 1:max_iter
global policy
global v

V_matrix = repeat(v, 1, N)
c_matrix = k_grid.^α .+ (1-δ).*k_grid .- transpose(k_grid);
v, policy_idx = find_maximiser(u.(c_matrix) + β*V_matrix);

policy_new = k_grid[policy_idx];

display(policy_new)
display(policy)
c_opt = k_grid.^α+(1-δ)*k_grid-policy;

Q = spzeros(N,N);
for i=1:N
Q[i, policy_idx[i]] = 1;
end

v_new = (Matrix(I,N,N)-β*Q) \ u.(c_opt);
diff = maximum(abs.(policy_new - policy));
println("Iteration \$iter, Error = \$diff")

v = v_new;
policy = policy_new;

if diff < tol
println("Solved at Iteration \$iter, Error = \$diff")
return policy, policy_idx
break
end

end
``````

It gives me the same answer if I wrap the whole code in `function foo() ... end` and delete the `global` declarations. Both in a function and in global scope it is converging at `Iteration 2, Error = 0.0` (not iteration 1), so I can’t replicate the result you claim.

PS. Just do `I - β*Q` — no need to construct a `Matrix`, which will discard your sparse data structure in favor of a dense matrix.

PPS. No need for the `repeat` call: just do `u.(c_matrix) .+ β.*v` and the broadcast property of `.+` should do the repetition implicitly.

PPPS. Looks like you want `global diff` in the loop too, since you write to the `diff` variable.

5 Likes

Ah sorry! What I mean by it works with functions is that I have a different implementation written with functions as follows. It’s not a “clone” of the code above but the basic structure is the same.

``````using Plots
using Parameters
using LinearAlgebra

u(c) =  c <= 0 ? -Inf : log(c);

Model = @with_kw (β = 0.95,
α = 0.65,
δ = 1,
N = 5,
k_ss = ((1/β-1+δ)/α)^(1/(α-1)),
cover = 0.1,
k_grid = range((1-cover)*k_ss, (1+cover)*k_ss, length=N),
T = 1000,
);

function find_maximiser(A)
v, B = findmax(A)
B = getindex(B,2)
return v, B
end

function solve(m, v; tol=1e-6, max_iter=500)
@unpack α, β, δ, N, k_grid = m
dr = zeros(N);
v1 = zeros(N)
crit = Inf
kp0 = k_grid
iter = 0
local kp
while crit>tol && iter <= max_iter

for i=1:N
c = k_grid[i].^α .+ (1-δ).*k_grid[i] .- k_grid;
v1[i], dr[i] = findmax(u.(c) + β*v);
end

kp = k_grid[Int.(dr)];
c = k_grid.^α+(1-δ)*k_grid-kp;

Q = zeros(N,N);
for i=1:N
Q[i, Int(dr[i])] = 1;
end

Tv = (Matrix(I,N,N) - β*Q)\u.(c);
crit= maximum(abs.(kp-kp0));
println("Iteration \$iter, Error = \$crit")
v = Tv;
kp0 = kp;
iter += 1
end
return kp
end

m = Model();
v = zeros(m.N)
pol = solve(m, v);

plt = plot(m.k_grid, pol, label="PFI")
plt = plot!(m.k_grid, m.k_grid, color =:black, linestyle =:dash, label="")
plt = plot!(xlabel="Capital Today", ylabel="Capital Tomorrow")

``````