diff --git a/HARK/ConsumptionSaving/ConsIlliquidAssetModel.py b/HARK/ConsumptionSaving/ConsIlliquidAssetModel.py index 478c044ef..06c94dc55 100644 --- a/HARK/ConsumptionSaving/ConsIlliquidAssetModel.py +++ b/HARK/ConsumptionSaving/ConsIlliquidAssetModel.py @@ -34,6 +34,22 @@ ############################################################################### +class ConFromAssetsFunction(MetricObject): + """ + A trivial class for wrapping an "assets function" in a consumption function. + """ + def __init__(self, aFunc): + self.aFunc = aFunc + + def __call__(self,x,y): + return x - self.aFunc(x,y) + + def derivativeX(self,x,y): + return 1. - self.aFunc.derivativeX(x,y) + + def derivativeY(self,x,y): + return -self.aFunc.derivativeY(x,y) + class DepositFunction(MetricObject): """ @@ -134,13 +150,16 @@ class SpecialMargValueFunction(MetricObject): None """ - distance_criteria = ["dFunc", "dvdk", "dvdb", "Penalty"] + distance_criteria = ["dFunc", "dvdk", "dvdb", "cFunc_Deposit", "cFunc_Withdraw", "Penalty", "CRRA"] - def __init__(self, dFunc, dvdk, dvdb, Penalty): + def __init__(self, dFunc, dvdk, dvdb, cFunc_Deposit, cFunc_Withdraw, Penalty, CRRA): self.dFunc = dFunc self.dvdk = dvdk self.dvdb = dvdb + self.cFunc_Deposit = cFunc_Deposit + self.cFunc_Withdraw = cFunc_Withdraw self.Penalty = Penalty + self.CRRA = CRRA def __call__(self, mNrm, nNrm): """ @@ -149,11 +168,21 @@ def __call__(self, mNrm, nNrm): """ dNrm = self.dFunc(mNrm, nNrm) withdraw = dNrm < 0.0 + deposit = dNrm > 0.0 dNrmAdj = dNrm * (1.0 + self.Penalty * withdraw) kNrm = mNrm - dNrm bNrm = np.maximum(nNrm + dNrmAdj, 0.0) dvdm = self.dvdk(kNrm, bNrm) dvdn = self.dvdb(kNrm, bNrm) + withdraw = np.logical_and(withdraw, bNrm > 0.) + if np.any(deposit): + cNrm_dep = self.cFunc_Deposit(bNrm[deposit]) + dvdm[deposit] = cNrm_dep**-self.CRRA + dvdn[deposit] = dvdm[deposit] + if np.any(withdraw): + cNrm_wdw = self.cFunc_Withdraw(bNrm[withdraw]) + dvdm[withdraw] = cNrm_wdw**-self.CRRA + dvdn[withdraw] = dvdm[withdraw] / (1.0 + self.Penalty) return dvdm, dvdn def dvdm(self, mNrm, nNrm): @@ -180,13 +209,15 @@ class IlliquidConsumerSolution(MetricObject): distance_criteria = ["cFunc", "dFunc"] - def __init__(self, MargValueFunc, cFunc, dFunc, mNrmMin, MPCmin=None, hNrm=None): + def __init__(self, MargValueFunc, cFunc, dFunc, mNrmMin, MPCmin=None, hNrm=None, cFunc_Deposit=None, cFunc_Withdraw=None): self.MargValueFunc = MargValueFunc self.cFunc = cFunc self.dFunc = dFunc self.mNrmMin = mNrmMin self.MPCmin = MPCmin self.hNrm = hNrm + self.cFunc_Deposit = cFunc_Deposit + self.cFunc_Withdraw = cFunc_Withdraw def make_basic_illiquid_solution_terminal(CRRA): @@ -535,19 +566,27 @@ def solve_one_period_basic_illiquid( # market resources functions dvdbNvrs = uFunc.derinv(EndOfPrd_dvdb, order=(1, 0)) cFunc_by_bNrm_list = [] + aFunc_by_bNrm_list = [] dvdbNvrsFunc_by_bNrm_list = [] for j in range(bCount): kNrm_temp = np.insert(kNrmNow[:, j] - BoroCnstNat[j], 0, 0.0) + aNrm_temp = np.insert(aNrmNow[:, j], 0, BoroCnstNat[j]) cNrm_temp = np.insert(cNrmNow[:, j], 0, 0.0) dvdbNvrs_temp = np.insert(dvdbNvrs[:, j], 0, 0.0) intercept_temp = MPCminNow * (hNrmNow + bNrmGrid[j] + BoroCnstNat[j]) cFunc_by_bNrm_list.append( LinearInterp(kNrm_temp, cNrm_temp, intercept_temp, MPCminNow) ) + aFunc_by_bNrm_list.append( + LinearInterp(kNrm_temp, aNrm_temp) + ) dvdbNvrsFunc_by_bNrm_list.append(LinearInterp(kNrm_temp, dvdbNvrs_temp)) - cFuncUncBase = LinearInterpOnInterp1D(cFunc_by_bNrm_list, bNrmGrid) - cFuncUnc = VariableLowerBoundFunc2D(cFuncUncBase, BoroCnstNatFunc) + #cFuncUncBase = LinearInterpOnInterp1D(cFunc_by_bNrm_list, bNrmGrid) + #cFuncUnc = VariableLowerBoundFunc2D(cFuncUncBase, BoroCnstNatFunc) cFuncCnstBase = IdentityFunction(i_dim=0, n_dims=2) + aFuncUncBase = LinearInterpOnInterp1D(aFunc_by_bNrm_list, bNrmGrid) + aFuncUnc = VariableLowerBoundFunc2D(aFuncUncBase, BoroCnstNatFunc) + cFuncUnc = ConFromAssetsFunction(aFuncUnc) if (BoroCnstArt is not None) and (BoroCnstArt > -np.inf): borrowing_constraint = LinearInterp( [0.0, 1.0], [BoroCnstArt, BoroCnstArt - 1 / (1.0 + IlqdPenalty)] @@ -660,8 +699,10 @@ def solve_one_period_basic_illiquid( plt.plot(bNrmGrid, OptZeroDeposit, "-b") plt.plot(bNrmGrid, OptZeroWithdraw, "-r") + plt.xlabel('illiquid assets bNrm') + plt.ylabel('liquid assets kNrm') plt.show() - + # Construct linear interpolations of the boundaries of the region of inaction InactionUpperBoundFunc = LinearInterp(bNrmGrid, OptZeroDeposit) InactionLowerBoundFunc = LinearInterp(bNrmGrid, OptZeroWithdraw) @@ -673,7 +714,19 @@ def solve_one_period_basic_illiquid( # Construct the "optimal kash on hand when withdrawing" function TotalWealthAdj = OptZeroWithdraw + bNrmGrid / (1.0 + IlqdPenalty) CashFunc_Withdraw = LinearInterp(TotalWealthAdj, OptZeroWithdraw) - + + # Construct simplified consumption functions when withdrawing or depositing + cNrm_dep = cFuncNow(OptZeroDeposit, bNrmGrid) + cNrm_wdw = cFuncNow(OptZeroWithdraw, bNrmGrid) + cFunc_Deposit = LinearInterp(bNrmGrid, cNrm_dep) + cFunc_Withdraw = LinearInterp(bNrmGrid, cNrm_wdw) + + plt.plot(bNrmGrid, cNrm_dep, "-b") + plt.plot(bNrmGrid, cNrm_wdw, "-r") + plt.xlabel('illiquid assets bNrm') + plt.ylabel('consumption cNrm') + plt.show() + # Construct the deposit function as a special representation dFuncNow = DepositFunction( InactionLowerBoundFunc, @@ -686,7 +739,7 @@ def solve_one_period_basic_illiquid( # Construct the marginal value function MargValueFuncNow = SpecialMargValueFunction( - dFuncNow, dvdkFunc, dvdbFunc, IlqdPenalty + dFuncNow, dvdkFunc, dvdbFunc, cFunc_Deposit, cFunc_Withdraw, IlqdPenalty, CRRA, ) # Package and return the solution object @@ -697,6 +750,8 @@ def solve_one_period_basic_illiquid( mNrmMin=mNrmMinNow, MPCmin=MPCminNow, hNrm=hNrmNow, + cFunc_Deposit=cFunc_Deposit, + cFunc_Withdraw=cFunc_Withdraw, ) return solution_now @@ -740,7 +795,7 @@ def solve_one_period_basic_illiquid( "bNrmMin": 0.0, # Minimum end-of-period "assets above minimum" value "bNrmMax": 16, # Maximum end-of-period "assets above minimum" value "bNrmNestFac": 1, # Exponential spacing for bNrmGrid - "bNrmCount": 300, # Number of points in the grid of "assets above minimum" + "bNrmCount": 96, # Number of points in the grid of "assets above minimum" "bNrmExtra": None, # Additional other values to add in grid (optional) } @@ -817,11 +872,18 @@ def cFuncSimple(self, t, mNrm, nNrm): cFunc = self.solution[t].cFunc dNrm = dFunc(mNrm, nNrm) kNrm = mNrm - dNrm - temp = 1.0 + self.IlqdPenalty * (dNrm < 0.0) + withdraw = dNrm < 0.0 + deposit = dNrm > 0.0 + temp = 1.0 + self.IlqdPenalty * withdraw bNrm = np.maximum(nNrm + temp * dNrm, 0.0) - # shift = cFunc.functions[0].lowerBound(bNrm) - # cNrm = cFunc.functions[0].func._secret(kNrm-shift,bNrm) cNrm = cFunc(kNrm, bNrm) + withdraw = np.logical_and(withdraw, bNrm > 0.) + if np.any(deposit): + cNrm_dep = self.solution[t].cFunc_Deposit(bNrm[deposit]) + cNrm[deposit] = cNrm_dep + if np.any(withdraw): + cNrm_wdw = self.solution[t].cFunc_Withdraw(bNrm[withdraw]) + cNrm[withdraw] = cNrm_wdw return cNrm def MPCfuncSimple(self, t, mNrm, nNrm): @@ -845,25 +907,15 @@ def MPCfuncSimple(self, t, mNrm, nNrm): if __name__ == "__main__": MyType = BasicIlliquidConsumerType() - MyType.cycles = 2 + MyType.cycles = 3 MyType.solve() - B = 8.0 + B = 6.0 f = lambda x: MyType.MPCfuncSimple(0, x, B * np.ones_like(x)) g = lambda x: MyType.cFuncSimple(0, x, B * np.ones_like(x)) h = lambda x: MyType.solution[0].dFunc(x, B * np.ones_like(x)) z = lambda x: (x + MyType.solution[0].hNrm + B) * MyType.solution[0].MPCmin c = lambda x: MyType.solution[0].cFunc(x, B * np.ones_like(x)) - plot_funcs([c, z], 0, 80) - - # m= 0.33446644398923286 n= 0.0390526406103266 - # m= 0.4051718438763465 n= 0.0390526406103266 - # m= 0.44712405210883144 n= 0.0390526406103266 - # m= 0.4844309585942469 n= 0.0390526406103266 - # m= 0.5231544276029979 n= 0.0390526406103266 - # m= 0.5703958661408141 n= 0.0390526406103266 - # m= 0.6620836699431031 n= 0.0390526406103266 - # m= -0.1852923544627234 n= 0.03462837348818072 - # m= 0.5091270511584698 n= 0.03462837348818072 - # m= 0.5464339576438852 n= 0.03462837348818072 + plot_funcs([g], -7.5, 12) + diff --git a/HARK/ConsumptionSaving/ConsIndShockModel.py b/HARK/ConsumptionSaving/ConsIndShockModel.py index 61b7725cc..3ed31892c 100644 --- a/HARK/ConsumptionSaving/ConsIndShockModel.py +++ b/HARK/ConsumptionSaving/ConsIndShockModel.py @@ -2083,19 +2083,19 @@ def make_euler_error_func(self, mMax=100, approx_inc_dstn=True): IncShkDstn = self.IncShkDstn[0] else: TranShkDstn = MeanOneLogNormal(sigma=self.TranShkStd[0]).discretize( - N=200, + N=100, method="equiprobable", - tail_N=50, + tail_N=25, tail_order=1.3, tail_bound=[0.05, 0.95], ) TranShkDstn = add_discrete_outcome_constant_mean( - TranShkDstn, self.UnempPrb, self.IncUnemp + TranShkDstn, p=self.UnempPrb, x=self.IncUnemp ) PermShkDstn = MeanOneLogNormal(sigma=self.PermShkStd[0]).discretize( - N=200, + N=100, method="equiprobable", - tail_N=50, + tail_N=25, tail_order=1.3, tail_bound=[0.05, 0.95], ) @@ -2118,12 +2118,12 @@ def make_euler_error_func(self, mMax=100, approx_inc_dstn=True): aNowGrid = mNowGrid - cNowGrid # Tile the grids for fast computation - ShkCount = IncShkDstn[0].size + ShkCount = IncShkDstn.pmv.size aCount = aNowGrid.size aNowGrid_tiled = np.tile(aNowGrid, (ShkCount, 1)) - PermShkVals_tiled = (np.tile(IncShkDstn[1], (aCount, 1))).transpose() - TranShkVals_tiled = (np.tile(IncShkDstn[2], (aCount, 1))).transpose() - ShkPrbs_tiled = (np.tile(IncShkDstn[0], (aCount, 1))).transpose() + PermShkVals_tiled = (np.tile(IncShkDstn.atoms[0], (aCount, 1))).transpose() + TranShkVals_tiled = (np.tile(IncShkDstn.atoms[1], (aCount, 1))).transpose() + ShkPrbs_tiled = (np.tile(IncShkDstn.pmv, (aCount, 1))).transpose() # Calculate marginal value next period for each gridpoint and each shock mNextArray = (