Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

One-asset HANK: Fail to solve steady state #10

Open
chicoutimy opened this issue Apr 7, 2022 · 4 comments
Open

One-asset HANK: Fail to solve steady state #10

chicoutimy opened this issue Apr 7, 2022 · 4 comments

Comments

@chicoutimy
Copy link

chicoutimy commented Apr 7, 2022

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_steady_state
_ = 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
Copy link

caromoree commented Jun 29, 2022

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
Copy link
Collaborator

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
Copy link

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!

mikkelpm added a commit to mikkelpm/stderr_calibration_python that referenced this issue Dec 31, 2022
Update hank.py - solution to steady state convergence issue proposed by @lukas-hack at shade-econ/sequence-jacobian#10
mikkelpm referenced this issue in mikkelpm/stderr_calibration_python Dec 31, 2022
@bbardoczy
Copy link
Collaborator

bbardoczy commented Jul 19, 2024

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?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

4 participants