**Example**
Consider the system $\begin{align} \frac{d}{dx} x_{1} &= f_{1} (x_{1},x_{2}) = x_{2}, \\ \frac{d}{dt } x_{2} &= f_{2}(x_{1},x_{2})= 3x_{1}^{2}-1. \end{align}$
# Stationary Points
Stationary points $(\bar{x}_{2}, \bar{x}_{1})$ have to satisfy $0=f_{1}(\bar{x}_{1},\bar{x}_{2})=\bar{x}_{2}$ and $0=f_{2}(\bar{x}_{1},\bar{x}_{2})= 3\bar{x}_{1}^{2}-1$ so the set of stationary points is given by $\left\{ \left( \sqrt{ \frac{1}{3} }, 0 \right) , \left( \sqrt{ -\frac{1}{3} } , 0 \right) \right\}$
# [[Direction Field]]
```run-python
import matplotlib.pyplot as plt
import numpy as np
from scipy.integrate import odeint
def f(x,t):
x1, x2 = x
f1 = x2
f2 = 3 * x1**2 - 1
return [f1,f2]
# generating a meshgrid, the points in which we will dra arrows
x1min, x1max, x1step = -1, 1, 0.1
x2min, x2max, x2step = -1, 1, 0.1
x1points = np.arange(x1min,x1max+x1step,x1step)
x2points = np.arange(x2min,x2max+x2step,x2step)
x1, x2 = np.meshgrid(x1points, x2points)
Df1, Df2 = f([x1,x2],0) #input for t doesn't matter since the equation is autonomous
fig, ax = plt.subplots(figsize=(9,9))
ax.quiver(x1,x2,Df1,Df2)
ax.axis([-2.0,2.0,-2.0,2.0])
ax.set_aspect('equal')
plt.show()
```
# Phase Portrait
```run-python
# time array for solution
tStart = 0.0
tEnd = 5.0
tInc = 0.05
t = np.arange(tStart,tEnd,tInc)
figPP, axPP = plt.subplots(figsize=(9,9))
# direction field
axPP.quiver(x1,x2,Df1,Df2)
# a couple of solutions
for p in np.arange(-0.5,0.75,0.25):
# compute numerical solution using forward-euler method
x0 = [p,0]
xsol = odeint(f,x0,t)
axPP.plot(xsol[:,0],xsol[:,1]) #,'.',ms=1)
# setting figure properties
axPP.set_xlabel('$x_1