Help. Some quantum stuff going on in my code. Probably related to scoping

I need help with code below. It is a value function iteration code–probably inefficient–but anyway. I have two vectors TV and v, both declared outside two loops: outer while loop and inner for loop. After the first pass of the while loop I assign v =TV following some computations using the for loop. TV here is the new approximation to the true value function, and v the old approximation. The problem is that on the second pass of the while loop, v and TV appear linked so that changing TV immediately changes v before I make the assignment when the for loop completes.

Probably easier to see with the code. To me this seems to be some sort of quantum stuff in my code, but it is probably related to scoping. Just can’t figure out how to correct it. I have added some display statements in the code below to illustrate where the problem lies.

start = time()

function ramsey()
    σ = 1.5;                # utility parameter: IIES
    δ = 0.1;                # depreciation rate
    β = 0.95;               # discount factor
    α = 0.30;               # capital elasticity in output
    
    nbk = 10;             # size of grid, i.e., number of data points in grid
    crit = 1.0;               # convergence criterion
    epsi = 1e-6;            # convergence parameter
    
    ks = ((1-β*(1-δ))/(α*β))^(1/(α-1));     # steady state capital stock
    
    dev = 0.9;              # maximum deviation from steady state
    kmin = (1-dev)*ks       # lower bound on grid
    kmax = (1+dev)*ks       # upper bound on grid
    
    dk = (kmax - kmin)/(nbk-1)  # grid increment
    kgrid = collect(kmin:dk:kmax);       # build grid
    
    v = zeros(nbk,1);
    TV = zeros(nbk,1);
    dr = zeros(Int,nbk,1);
    util = zeros(nbk,1);
    
    while crit > epsi
        for i in 1:nbk
            #for i in eachindex(kgrid)
            #i==1 && display([v,TV]);
            #if i==1; display([v,TV]) end;
            #i==1 ? display([v,TV]) : "";
            
            # compute indices for each positive consumption
            #tmp = (kgrid[i]^α + (1-δ)*kgrid[i] - kmin);
            tmp = (kgrid[i]^α + (1-δ)*kgrid[i] - kmin);
            
            #imax = min(floor(Int,tmp/dk)+1,nbk);
            imax = min(floor(Int,tmp/dk)+1,nbk);
            
            # consumption and utility
            #c = kgrid[i]^α + (1-δ)*kgrid[i] .- kgrid[1:imax];
            c = kgrid[i]^α + (1-δ)*kgrid[i] .- kgrid[1:imax];
            
            #util = (c.^(1-σ) .-1) ./(1-σ);
            util = (c.^(1-σ) .-1) ./(1-σ);
            
            # find value function
            #(TV[i], dr[i]) = findmax(util +β*v[1:imax]);
            #(i == 1 || i == 2) && display([i, "Sam1", v,TV]);
            (i   == 1) && display([i, "Sam1", v,TV]);
            (TV[i], dr[i]) = findmax(util + β.*v[1:imax]);
            #(i == 1 || i == 2) && display([i, "Sam2", v,TV]);
            (i == 1) && display([i, "Sam2", v,TV]);
            #display(v)
        end
        ΔV = TV - v;
        v = TV;
        crit = maximum(abs.(ΔV));
        display(["Sam3", v,TV])
        #TV .= 1            # update value function
    end;
    
    
    # final solution
    kp = kgrid[dr];
    c = kgrid.^α + (1-δ).*kgrid - kp;
    util = (c.^(1-σ) .- 1) ./(1-σ);
    return kgrid, kp, c, util;
end


#kgrid, kp, c, util = ramsey()
kgrid, kp, c, util = ramsey();
elapsed = time() - start

#using Plots;
#plot(kgrid,c)

This is expected. The assignment v = TV does not copy, but makes both names refer to the same data.
You can use v .= TV to copy the conents of vector TV into vector v. This is different from e.g. v = copy(TV) or v = Array(TV) which creates a new vector.

2 Likes

Thanks Andreas, now I KNOW!! Works perfectly and now I understand another difference between “.=” and “=”

You’re welcome! The technical term is “variable binding” which you can search for. (I’m afraid I could not find a simple explanation in the manual)

1 Like