ODEs

We will use this notebook to define the differential equations that will be solved.

Each function representing a differential equation or set of differential equations should take a state vector and time array as its first two arguments, all other arguments relating to the system should follow after.

The function below describes the equation $\frac{dy}{dt}=y$ which has the analytical solution $y(t)=y_0e^t$


In [1]:
def Exponential(y, t, args=None):
    dydt=y
    return(dydt)

The function below describes the equation $\dfrac{dy}{dt}=\dfrac{-y}{\tau}$ which has the analytical solution $y(t)=y_0e^{-\dfrac{t}{\tau}}$


In [2]:
def RadioactiveDecay(y, t, tau):
    dydt=-1.0*y/tau
    return(dydt)

The function below describes the equation $\dfrac{dy}{dt}=t-y^2$


In [3]:
def Nonlinear1(y,t,args=None):
    dydt=t-np.square(y)
    return(dydt)

The function below describes the equation $\dfrac{dy}{dt}=-y^3+\sin{t}$


In [4]:
def Nonlinear2(y,t,args=None):
    dydt=-y**3 + np.sin(t)
    return(dydt)

The function below represents the motion of a simple damped pendulum which is described by the following equation. $$ \ddot\theta(t) + \alpha \dot\theta(t) + \frac{g}{L} \sin(\theta(t))=0$$ This is a second order system and must be split into linear parts to calculate a numerical solution. This splitting works as follows: $$ \begin{align*} \dot \theta(t) &= \omega(t) \\ \dot \omega(t) &= -\alpha \omega(t) -\beta \sin(\theta(t)) \end{align*} $$

Where $\beta=\dfrac{g}{L}$

In this case the state vector $y$ consists of the angular position and angular velocity of the pendulum bob. $$y = [\theta, \omega]$$


In [5]:
def Pendulum(y,t,alpha, beta):
    theta, omega = y
    dydt = np.array([omega, -alpha*omega - beta*np.sin(theta)])
    return dydt

The function below describes a simple pendulum which has had a small angle approximation in normalised units

$$\dfrac{d^2x}{dt^2} + x = 0$$

This can be split as follows and solved numerically:

$$ \begin{align*} \dfrac{dx}{dt}&=y\\ \dfrac{dy}{dt}&=-x \end{align*} $$

In [6]:
def SimplePendulum(y,t,b=0,omega=1):
    x,v  = y
    dydt = np.array([v,-b*v-(omega**2)*x])
    return dydt

The cell below defines the differential form of the Logistic equation:

$\dfrac{dx}{dt}=rx(1-x)$


In [7]:
def LogisticEquation(y,t,r,k=1):
    dydt=r*y*(1-y/k)
    return(dydt)

The cell below defines the Lotka-Volterra Model which describes population dynamics

$$ \begin{align*} \dfrac{dx}{dt}&=ax-bxy\\ \dfrac{dy}{dt}&=cxy-dy \end{align*} $$

In [8]:
def LotkaVoltera(f,t,a,b,c,d):
    x,y=f
    dxdt=a*x-b*x*y
    dydt=c*x*y-d*y
    dfdt=np.array([dxdt,dydt])
    return(dfdt)

In [9]:
def Nonlinear3(f,t):
    x,y=f
    dxdt=y
    dydt=-x+(1-x**2)*y
    dfdt=np.array([dxdt,dydt])
    return(dfdt)

In [10]:
def NonlinearSin(x,t):
    return(np.sin(x))

In [11]:
def Budworm(x,t,r,k):
    dxdt=r*x*(1-x/k)-(x**2/(1+x**2))
    return(dxdt)

In [12]:
def Nonlinear4(f,t):
    x,y=f
    dxdt=x-y
    dydt=(x*x)-4
    dfdt=np.array([dxdt,dydt])
    return(dfdt)

In [13]:
def Nonlinear5(f,t):
    x,y=f
    dxdt=-x-y
    dydt=(x*x)-4
    dfdt=np.array([dxdt,dydt])
    return(dfdt)

In [ ]: