This is a first go at using a computer algebra system (SymPy) to explore the formal steady state of a Stock-Flow Consistent model. For this, I use a simplified five-equation version of Godley and Lavoie’s (2007) Model SIM. For more information, see my previous posts here and here.


Assume a closed economy with no investment such that national income, , in time consists of household final consumption expenditure, , and government expenditure, :

\begin{equation} \tag{1} Y_t = C_t + G_t \end{equation}

Government expenditure is exogenous. Consumption depends on current disposable income, , and accumulated wealth, (from the previous period):

\begin{equation} \tag{2} C_t = \alpha_1 \cdot YD_t + \alpha_2 \cdot H_{t-1} \end{equation}

and are the marginal propensities to consume out of disposable income and wealth , respectively.

Disposable income is national income minus taxes, , which are levied at a fixed proportion, :

\begin{equation} \tag{3} YD_t = Y_t - T_t \end{equation}

\begin{equation} \tag{4} T_t = \theta \cdot Y_t \end{equation}

Households accumulate savings from the difference between their disposable income and their expenditure :

\begin{equation} \tag{5} H_t = H_{t-1} + YD_t - C_t \end{equation}

Numerical solution

We can use FSIC to define and then solve the model numerically (see my earlier post for details of the syntax):

import fsic

# Define the model's economic logic in a string
script = '''
Y = C + G
C = {alpha_1} * YD + {alpha_2} * H[-1]
YD = Y - T
T = {theta} * Y
H = H[-1] + YD - C

# Translate into a set of symbols and a model class definition in one go
SIM = fsic.build_model(fsic.parse_model(script))

Simulate the model over a long enough period for it to achieve its stationary state (at which the values of the variables in levels are constant):

model = SIM(range(1, 100 + 1))  # Solve from Period 1 to 100

# Set propensities to consume
model.alpha_1 = 0.6
model.alpha_2 = 0.4

# Fiscal policy
model.G[1:] = 20   # Government spending
model.theta = 0.2  # Rate of income tax


We can then reproduce the main elements of Table 3.4 from Godley and Lavoie (2007) (‘The impact of $20 of government expenditures, with perfect foresight’).

While I’ve made no changes to FSIC for this post, I have added a new module to accompany the core library: fsictools. This provides the to_dataframe() function from an earlier post.

from fsictools import to_dataframe

# Convert the contents of the model to a `pandas` DataFrame
results = to_dataframe(model)

# Extract selected variables of interest, adding a new column for the change in
# wealth
extract = results[['G', 'Y', 'T', 'YD', 'C', 'H']].copy()
extract['D(H)'] = extract['H'].diff()

# Transpose, and select periods, to reproduce the main elements of Table 3.4
table_3_4 = extract.T
table_3_4[[1, 2, 3, 100]].round(1)
$$1$$ $$2$$ $$3$$ $$\infty$$
$$G$$ 0.0 20.0 20.0 20.0
$$Y$$ 0.0 38.5 47.9 100.0
$$T$$ 0.0 7.7 9.6 20.0
$$YD$$ 0.0 30.8 38.3 80.0
$$C$$ 0.0 18.5 27.9 80.0
$$H$$ 0.0 12.3 22.7 80.0
$$\Delta H$$ 0.0 12.3 10.4 0.0

By the one-hundredth period, the model has achieved a stationary state at which the values of the variables in levels are constant. The period-on-period change in household wealth is zero.

Analytical solution

Model SIM is small and thus amenable to more formal algebraic analysis of its properties. The remainder of this post replicates key results from Section 3.5 of Godley and Lavoie (2007) using SymPy. I’m new to computer algebra systems so this is both a first go and a note of my early learnings.

Define the system of equations in SymPy. Below, I use 1 as a suffix to indicate the one-period lag on wealth .

from sympy import Symbol, Eq        # Term and equation objects
from sympy import symbols, sympify  # Helpers to generate SymPy objects
from sympy import factor, linsolve  # Expression manipulation and system solution

# Define symbols for the endogenous variables
Y, C, YD, T, H = symbols('Y, C, YD, T, H')

# Define the individual equations with left-hand and right-hand
# side expressions
output            = Eq(lhs=Y,  rhs=sympify('C + G'))
consumption       = Eq(lhs=C,  rhs=sympify('alpha_1 * YD + alpha_2 * H1'))
disposable_income = Eq(lhs=YD, rhs=sympify('Y - T'))
taxes             = Eq(lhs=T,  rhs=sympify('theta * Y'))
wealth            = Eq(lhs=H,  rhs=sympify('H1 + YD - C'))

# Assemble the system of equations as a dictionary, keyed by LHS
# variable (symbol)
system = {
    Y: output,
    C: consumption,
    YD: disposable_income,
    T: taxes,
    H: wealth,
# Display the system
for v in system.values():

\begin{equation} Y = C + G \end{equation}

\begin{equation} C = H_{1} \alpha_{2} + YD \alpha_{1} \end{equation}

\begin{equation} YD = - T + Y \end{equation}

\begin{equation} T = Y \theta \end{equation}

\begin{equation} H = - C + H_{1} + YD \end{equation}

At the model’s stationary state, the values of the variables are constant in levels. As in the above numerical solution, the change in wealth should be zero:

\begin{equation} \Delta H_t = H_t - H_{t-1} = 0 \end{equation}

\begin{equation} H_t = H_{t-1} = H_{t-2} = H_{t-3} = \dots \end{equation}

We can impose this on the system by replacing instances of H1 (the one-period lag of wealth) with H (the contemporaneous value of wealth):

stationary_state = {
    k: v.subs(Symbol('H1'), Symbol('H'))
    for k, v in system.items()
# Display the system
for v in stationary_state.values():

\begin{equation} Y = C + G \end{equation}

\begin{equation} C = H \alpha_{2} + YD \alpha_{1} \end{equation}

\begin{equation} YD = - T + Y \end{equation}

\begin{equation} T = Y \theta \end{equation}

\begin{equation} H = - C + H + YD \end{equation}

We can solve for the stationary state with linsolve() and inspect the values:

values = list(linsolve(list(stationary_state.values()),    # System of equations
                       list(stationary_state.keys())))[0]  # Unknowns to solve for
solution = dict(zip(stationary_state.keys(), values))
# Display the solution at its stationary state
for k, v in solution.items():
    display(Eq(k, factor(v)))

\begin{equation} Y = \frac{G}{\theta} \end{equation}

\begin{equation} C = - \frac{G \left(\theta - 1\right)}{\theta} \end{equation}

\begin{equation} YD = - \frac{G \left(\theta - 1\right)}{\theta} \end{equation}

\begin{equation} T = G \end{equation}

\begin{equation} H = \frac{G \left(\alpha_{1} - 1\right) \left(\theta - 1\right)}{\alpha_{2} \theta} \end{equation}

As in Equation 3.15 of Godley and Lavoie (2007), the stationary state for national income (the first result, ) is:

\begin{equation} Y^{\star} = \frac{G}{\theta} \end{equation}

represents the fiscal stance: the ratio of government expenditure to its income share (also as in Godley and Cripps, 1983).

The expressions for household final consumption expenditure and disposable income are the same. The two are equal. This is consistent with (but, here, also follows from) there being no change in wealth at the stationary state. With slight rearrangement:

\begin{equation} C^{\star} = YD^{\star} = G \cdot \frac{1 - \theta}{\theta} \end{equation}

At the stationary state, tax revenues match government expenditures such that government debt (the mirror image of household wealth) is constant.

Finally, the level of household wealth (or, conversely, government debt) is, again with some rearrangement:

\begin{equation} H^{\star} = \frac{1 - \alpha_1}{\alpha_2} \cdot G \cdot \frac{1 - \theta}{\theta} = \alpha_3 \cdot G \cdot \frac{1 - \theta}{\theta} = \alpha_3 \cdot YD^{\star} \end{equation}

\begin{equation} \alpha_3 = \frac{1 - \alpha_1}{\alpha_2} \end{equation}

Update 05/01/2020

Having solved for the endogenous variables at the model’s stationary state, we can substitute in the numerical values of the exogenous variables to check that the analytical solution matches our numerical solution above.

# Set values for the exogenous variables using a dictionary:
#  - keys: variables as `Symbol` objects
#  - values: numerical values
substitutions = {Symbol(k): v
                 for k, v in {'alpha_1': 0.6, 'alpha_2': 0.4,
                              'G': 20.0, 'theta': 0.2}.items()}

# Display the solution's values at its stationary state for the given inputs
for k, v in solution.items():
    display(Eq(k, v.subs(substitutions)))

\begin{equation} \displaystyle Y = 100.0 \end{equation}

\begin{equation} \displaystyle C = 80.0 \end{equation}

\begin{equation} \displaystyle YD = 80.0 \end{equation}

\begin{equation} \displaystyle T = 20.0 \end{equation}

\begin{equation} \displaystyle H = 80.0 \end{equation}

See here for this post as a Jupyter notebook along with supporting Python code.


Godley, W., Cripps, F. (1983) Macroeconomics, Oxford University Press

Godley, W., Lavoie, M. (2007) Monetary economics: An integrated approach to credit, money, income, production and wealth, Palgrave Macmillan