One-asset HANK: Fail to solve steady state

[chicoutimy]

High, I am following Tutorial 3 to solve the one-asset HANK model. The code runs smoothly until line 8, but when I go ahead and run line 9

calibration = {'eis': 0.5, 'frisch': 0.5, 'rho_e': 0.966, 'sd_e': 0.5, 'nE': 7, 'amin': 0.0, 'amax': 150, 'nA': 500, 'Y': 1.0, 'Z': 1.0, 'pi': 0.0, 'mu': 1.2, 'kappa': 0.1, 'rstar': 0.005, 'phi': 1.5, 'B': 5.6}

unknowns_ss = {'beta': 0.986, 'vphi': 0.8} targets_ss = {'asset_mkt': 0, 'labor_mkt': 0}

ss0 = hank_ss.solve_steady_state(calibration, unknowns_ss, targets_ss, solver="hybr")

I get the following error message:

Traceback (most recent call last): File "C:\Users\user\anaconda3\lib\site-packages\IPython\core\interactiveshell.py", line 3369, in run_code exec(code_obj, self.user_global_ns, self.user_ns) File "", line 8, in <cell line: 8> ss0 = hank_ss.solve_steady_state(calibration, unknowns_ss, targets_ss, solver="hybr") File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\block.py", line 159, in solve_steadystate = solve_for_unknowns(residual, unknowns, solver, opts['solver_kwargs'], File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\support\steady_state.py", line 185, in solve_for_unknowns result = opt.root(residual_f, initial_values, File "C:\Users\user\anaconda3\lib\site-packages\scipy\optimize_root.py", line 187, in root sol = _root_hybr(fun, x0, args=args, jac=jac, *options) File "C:\Users\user\anaconda3\lib\site-packages\scipy\optimize\minpack.py", line 226, in _root_hybr shape, dtype = _check_func('fsolve', 'func', func, x0, args, n, (n,)) File "C:\Users\user\anaconda3\lib\site-packages\scipy\optimize\minpack.py", line 24, in _check_func res = atleast_1d(thefunc(((x0[:numinputs],) + args))) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\block.py", line 156, in residual ss.update(self.steady_state(ss, dissolve=dissolve, options=options, kwargs)) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\block.py", line 47, in steady_state return self.M @ self._steady_state(self.M.inv @ calibration, dissolve=dissolve, File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\combined_block.py", line 59, in _steady_state outputs = block.steady_state(ss, dissolve=inner_dissolve, kwargs) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\block.py", line 50, in steady_state return self.M @ self._steady_state(self.M.inv @ calibration, own_options) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\het_block.py", line 94, in _steady_state ss = self.backward_steady_state(ss, tol=backward_tol, maxit=backward_maxit) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\blocks\het_block.py", line 190, in backward_steady_state ss.update(self.backward_fun(ss)) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\utilities\function.py", line 68, in call return self.outputs.dict_from(make_tuple(self.f(input_dict))) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\hetblocks\hh_labor.py", line 33, in hh c[iconst], n[iconst] = solve_cn(we[iconst[0]], File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\hetblocks\hh_labor.py", line 51, in solve_cn uc = solve_uc(w, T, eis, frisch, vphi, uc_seed) File "C:\Users\user\anaconda3\lib\site-packages\sequence_jacobian\hetblocks\hh_labor.py", line 69, in solve_uc raise ValueError("Cannot solve constrained household's problem: No convergence after 30 iterations!") ValueError: Cannot solve constrained household's problem: No convergence after 30 iterations!

Do you know what is going on? Am I am doing something wrong? Thanks for your help and congrats on the wonderful toolbox! chicoutimy

[caromoree]

Hi! The code runs correctly in python 3.8 versions but i'm having the same problem as @chicoutimy with a 3.9 version. Still don't know why...

[bbardoczy]

Hi, Yes, it seems like numba.vectorize that we use when solving for the constrained households' labor supply decision does not work properly for some (new) versions of numba. It seems that this problem only arises on Windows machines (correct me if you have Mac or Linux).

There are two workarounds.

1. Use Python 3.8 with associated packages. If you need Python 3.9 for some other work, consider keeping multiple environments. Switching between these is basically costless using conda. See this link on how to do this.
2. Replace the vectorized implementation with an njit-ed loop over the gridpoints where the borrowing constraint binds. This will achieve similar speed at the cost of less elegant code.

Hope this helps!

[lukas-hack]

Hi,

I tried to look into this problem.

Background: I used Windows 11 and Python 3.8.15 as I downgraded following Bence's comment. It turned out that this did not resolve the problem. My old laptop with Python 3.7.3 did not have the problem. The update with which the problem came must be in between I suppose.

Findings:

• The vectorize decorator works fine. Replacing it with native python code (an extra outer loop) leaves the problem unchanged.

• It turns out that the problem is the following: the njit-ed netexp() function calls the also njit-ed cn() function. Somehow this leads to not well-defined behavior. I am not sure why but maybe a numba expert could tell. It is quite easy to see this problem in debug mode by manually going through the loop in solve_uc().

Solution: Just use the following netexp() function that does not call cn():

@njit
def netexp(log_uc, w, T, eis, frisch, vphi):
    uc = np.exp(log_uc)
    c = uc ** (-eis)
    n = (w * uc / vphi) ** frisch
    ne = c - w * n - T

    # c and n have elasticities of -eis and frisch wrt log u'(c)
    c_loguc = -eis * c
    n_loguc = frisch * n
    netexp_loguc = c_loguc - w * n_loguc
    return ne, netexp_loguc

The only (trivial) change I made is to copy and paste the first 3 lines instead of calling the cn() function.

Hope this helps!

[bbardoczy]

Hmm, interesting. I nest njit-ed functions all the time. But then again, I haven't experienced this issue in the first place.

Has @lukas-hack 's solution worked for everybody here since 2022?