1"""Full static block for one period: tasks + structures + land + energy + compute utilization + housing.23Given predetermined stocks (K_N, K_Q, K_E, K_H), the period's exogenous state (alpha_t, Z_t, eps_t) and4parameters, `static_period` returns all quantities and prices of the period except those that depend on5household consumption (housing rent), which `housing_rent` adds once total consumption is known.6`income_identity` checks the national-accounting identity GDP = factor payments.7"""8from __future__ import annotations910import numpy as np11from scipy.optimize import brentq1213from .tasks import Technology, technology_from_params141516def energy_supply(p, KE, gf=1.0):17 """Energy output; the resource endowment is measured in efficiency units R^E * gf."""18 return p.ZE * KE ** (1.0 - p.eta) * (p.RE * gf) ** p.eta192021def static_period(p, tech: Technology, KN, KQ, KE, KH, alpha_t, Z_t, eps_t, x0=None, gf=1.0):22 """gf = (1+g_Z)^t: Harrod-neutral efficiency factor applied to labor, land and the energy resource."""23 E = energy_supply(p, KE, gf) if p.nu_E > 0 else 0.024 H_eff = np.asarray(p.H, dtype=float) * gf25 T_eff = p.T_land * gf26 TH_eff = p.TH_land * gf2728 def run(u):29 Q = u * KQ30 EY = E - eps_t * Q if p.nu_E > 0 else 0.031 sol = tech.solve(Q, H_eff, KN, T_eff, EY, Z_t, alpha_t, x0=x0)32 return sol, EY3334 u = 1.035 if p.nu_E > 0 and KQ > 0:36 # feasibility: energy left for final production must be positive37 u_max = min(1.0, (1 - 1e-10) * E / (eps_t * KQ)) if eps_t * KQ > 0 else 1.038 sol, EY = run(u_max)39 net = sol["A"] - eps_t * sol["PE"]40 if net < 0.0:41 # capacity is idle until the marginal unit just covers its energy cost: bisection in log u42 lo = np.log(min(1e-12, 0.5 * u_max))43 hi = np.log(u_max)44 f_lo = (lambda s: s[0]["A"] - eps_t * s[0]["PE"])(run(np.exp(lo)))45 if f_lo <= 0:46 u = np.exp(lo)47 else:48 for _ in range(60):49 mid = 0.5 * (lo + hi)50 fm = (lambda s: s[0]["A"] - eps_t * s[0]["PE"])(run(np.exp(mid)))51 if fm > 0:52 lo = mid53 else:54 hi = mid55 if hi - lo < 1e-11:56 break57 u = np.exp(lo)58 sol, EY = run(u)59 else:60 u = u_max61 else:62 sol, EY = run(1.0)6364 Q = u * KQ65 PE = sol["PE"] if p.nu_E > 0 else 0.066 RQ = u * (sol["A"] - eps_t * PE) # net rental per unit of installed capacity67 RE = (1.0 - p.eta) * PE * E / KE if (p.nu_E > 0 and KE > 0) else 0.068 rentR = p.eta * PE * E if p.nu_E > 0 else 0.069 HS = KH ** (1.0 - p.nu_H) * TH_eff ** p.nu_H70 out = dict(Y=sol["Y"], w=sol["w"], A=sol["A"], RN=sol["RN"], RT=sol["RT"], PE=PE, RQ=RQ, RE=RE, rentR=rentR,71 util=u, Qtot=Q, E=E, EY=EY, EQ=eps_t * Q, HS=HS, L=sol["L"], Qi=sol["Q"], X=sol["X"], F=sol["F"],72 XB=sol["XB"], x=sol["x"], success=sol["success"], KN=KN, KQ=KQ, KE=KE, KH=KH, gf=gf,73 T_eff=T_eff, TH_eff=TH_eff, H_eff=H_eff,74 wage_bill=float(np.sum(sol["w"] * H_eff)))75 return out767778def housing_rent(p, per, C_total):79 """Rent per unit of housing services when total goods consumption is C_total (Cobb-Douglas share)."""80 RH = p.theta_H / (1.0 - p.theta_H) * C_total / per["HS"]81 per["RH"] = RH82 per["RHK"] = (1.0 - p.nu_H) * RH * per["HS"] / per["KH"] # rental of housing structures83 per["RTH"] = p.nu_H * RH * per["HS"] / per["TH_eff"] # residential land rent per efficiency unit84 per["GDP"] = per["Y"] + RH * per["HS"]85 return per868788def income_identity(p, per, markup_profit=0.0):89 """Residual of GDP - (wages + rentals + rents + profits). Zero in the model; must fail if accounting is broken."""90 inc = (per["wage_bill"] + per["RN"] * per["KN"] + per["RT"] * per["T_eff"] + per["RQ"] * per["KQ"]91 + per["RE"] * per["KE"] + per["rentR"] + per["RHK"] * per["KH"] + per["RTH"] * per["TH_eff"])92 # RQ already includes the markup profit when markup_Q > 1 (owners receive A - eps PE per unit)93 return per["GDP"] - inc949596def value_added_by_sector(p, per):97 """Sectoral value added (final goods, compute, energy, housing); sums to GDP."""98 final = per["Y"] - per["A"] * per["Qtot"] - per["PE"] * per["EY"]99 compute = per["A"] * per["Qtot"] - per["PE"] * per["EQ"]100 energy = per["PE"] * per["E"]101 housing = per["RH"] * per["HS"]102 return dict(final=final, compute=compute, energy=energy, housing=housing, total=final + compute + energy + housing)103