Geometric way of interpreting equations

$\frac{d}{dt}x = \sin x$

$t = \ln \frac{\csc x_0 + \cot x_0}{\csc x + \cot x}$

  • After a while, we end up with an implicit solution: $t = \ln \frac{\csc x_0 + \cot x_0}{\csc x + \cot x}$
  • Can we make out anything meaningful from this equation?
  • Can we extract some information from this equation without solving it?

Implicit plot

We can recast the equations as $\exp{t}\times (\csc x + \cot x) = (\csc x_0 + \cot x_0), $ and then look at the isoline or the contour of the equations.

  • By plotting the isoline in the t-x space, we can figure out the evolution of the equation in time implicitly
In [1]:
import numpy as np;
import matplotlib.pyplot as plt;
plt.rcParams.update({"text.usetex":True});
%config InlineBackend.figure_format = "svg"
from ipywidgets import interactive
In [2]:
def show_sol(x0=4.5):
    t = np.linspace(1e-6, 10, 100);
    x = np.linspace(1e-6, 2*np.pi, 100);
    [X, T] = np.meshgrid(x, t);
    F = np.exp(T)*(1/np.sin(X) + 1/np.tan(X));

    isovalue = (1/np.sin(x0) + 1/np.tan(x0))

    plt.contour(T, X, F, [isovalue]);
    plt.xlabel("t");
    plt.ylabel("X");
    plt.contour(T, X, X, [np.pi/2, 3*np.pi/2], colors='k')
    
w = interactive(show_sol, x0=(0.1, 6.2, 0.1));
In [3]:
w
In [4]:
x = np.linspace(0, 2*np.pi);
dxdt = np.sin(x)
plt.plot(x, dxdt)
plt.xlabel("x"); plt.ylabel("dx/dt")
plt.plot(x, 0*x);

We can solve the initial value problem using the function solve_ivp which can be found in scipy.integrate. A quick demonstration of the same is shown below.

In [5]:
from scipy.integrate import solve_ivp

def Ndot(t, N, K, r):
    return r*N*(1-N/K)

K = 3;
r = 0.3;
N0 = K/3;
t = 2;

sol = solve_ivp(Ndot, [0, 20], [N0], args=[K,r], rtol=1e-6);
plt.plot(sol.t, sol.y.T);
plt.plot(sol.t, K*np.ones(np.shape(sol.t)), '--k')
plt.xlabel("$t$"); plt.ylabel("N(t)");

We can vary more parameters to see their influence on the time evolution of the population.

In [6]:
def plot_logistic(K=4, r=0.3, N0 = 0.2, t_end =30):

    def Ndot(t, N, K, r):
        return r*N*(1-N/K)

    #K = 3;
    #r = 0.3;
    #N0 = K/3;
    
    sol = solve_ivp(Ndot, [0, t_end], [N0], args=(K,r), rtol=1e-6);
    plt.plot(sol.t, sol.y.T);
    plt.plot(sol.t, K*np.ones(np.shape(sol.t)), '--k'); 
    plt.plot(sol.t, K/2*np.ones(np.shape(sol.t)), '--r'); 
    plt.xlabel("$t$"); plt.ylabel("N(t)");
    
w = interactive(plot_logistic, K = (1, 10, 0.5), r = (0.1, 0.95, 0.1), N0 = (0.1, 20, 0.5), t_end=(10, 50, 10));
w
In [ ]: