Symbolic Case

You know the formula explicitly

$$f(x) = ... $$

In this case there are two options

  1. The integral can be solved analytically $\int f(x) dx = \text{Something you can write down}$
  2. The integral cannot be solved analytically

We can try sympy at first to see if the integral can be solved analytically

Part 1: "Solvable" Integrals

Example: Find $\int \sin^3(x) e^{-5x} dx$

Example: Find $\int \cos(bx)e^{-ax} dx$

Example: Find $\int \frac{(1+\sqrt{x})^{1/3}}{\sqrt{x}} dx$

Example (Definite) Find $\int_{0}^{\ln(4)}\frac{e^x}{\sqrt{e^{2x}+9}} dx$

Example (Improper) Find $ \int_{0}^{\infty} \frac{16 \tan^{-1}(x)}{1+x^2} dx $

Part 2: "Unsolvable" Integrals

Example: Find $\int_{1}^{2} e^{-\sin(x)} dx$

Won't run

So we use scipy's quad function to integrate numerically

Example: Find $\int_{0}^{2\pi} \frac{1}{(a-\cos(x))^2 + (b-\sin(x))^2} dx$

This won't run as this does not have an analytical solution. We must solve it numerically, and thus we need choose values of $a$ and $b$ each time

We can solve this for many different values of $a$ and $b$ quickly as follows:

Numerical Case

This is the case where one collects actual data and wants to find the integral

Example: Electric Pulses in particle physics experiments

Suppose this is current measured $I(t)$

Now in some detector experiments, energy is proportional to integrated current so $\text{Energy} \propto \int I(t) dt$

Example: Covid Data

Plot cases per day

Can find the cumulative number of cases by taking the integral

And plot

Monte Carlo Integration

First Order Ordinary Differential Equation

Air friction while falling

$$ \frac{dv}{dt} - \alpha v^2 + \beta = 0 \hspace{10mm} v(0) = 0$$

The first thing we need to do is write it in the form

$$\frac{dv}{dt} = f(t,v)$$

In other words, "derivative of v equals something that depends v and time". This is easy in this example:

$$\frac{dv}{dt} = \alpha v^2 - \beta$$

and note that there is no dependence on time in this particular example. We need to write this differential equation in python form

Solving differential equation

There are two main solvers in scipy

Method 1

Method 2

Look at the solution. It will become obvious why it is returned in this form once we deal with systems of ODEs (what these solvers are really meant for)

Coupled first order ODEs

$$ y_1' = y_1 + y_2^2 + 3x \hspace{10mm} y_1(0)=0$$$$ y_2' = 3y_1 + y_2^3 - \cos(x) \hspace{10mm} y_2(0)=0$$

Letting $S=(y_1, y_2)$ we need to write a function that returns $dS/dx = (dy_1/dx, dy_2/dx)$. The function $dS/dx$ can take in $S=(y_1, y_2)$ and $x$. This is like before, but in vector format

$$ \vec{S} = \begin{bmatrix} y_1 \\ y_2 \end{bmatrix} \hspace{10mm} \implies \hspace{10mm} \frac{d\vec{S}}{dx} = \vec{f}(x, \vec{S}) = \vec{f}(x, y_1, y_2) = \begin{bmatrix} y_1' \\ y_2' \end{bmatrix} = \begin{bmatrix} y_1 + y_2^2 + 3x\\ 3y_1 + y_2^3 - \cos(x) \end{bmatrix}$$

Second Order Differential Equation

Python does not have functions to directly solve second order ODEs.

Consider

$$\ddot{x} = -\dot{x}^2 + \sin(x)$$

We can convert this into two first order ODEs as follows:

Our two equations:

$$\dot{x} = v$$$$\dot{v} = -v^2 + \sin(x)$$

These are two coupled first order equations. They require an initial condition ($x_0$ and $v_0$)

Two coupled third order equations

$$\dddot{x_2}= -\ddot{x_1}^3 + \dot{x_2} + x_1 + \sin(t)$$

$$\dddot{x_1}= -2\dot{x_2}^2 + x_2$$

Define

$v_1 = \dot{x_1}$

$v_2 = \dot{x_2}$

$a_1 = \ddot{x_1} = \dot{v_1}$

These make up four differential equations. Then noting that $\dot{a_1} = \dddot{x_1}$ and $\dot{a_2} = \dddot{x_2}$ we get

$\dot{a_2} = -a_1^3 + v_2 + x_1 + \sin(t)$

$\dot{a_1} = -2v_2^2 + x_2$

Then

$$\vec{S} = \begin{bmatrix} x_1\\ v_1 \\ a_1 \\ x_2 \\ v_2 \\ a_2 \end{bmatrix} \hspace{10mm} \implies \hspace{10mm} \frac{d\vec{S}}{dt} = \begin{bmatrix} \dot{x_1}\\ \dot{v_1} \\ \dot{a_1} \\ \dot{x_2} \\ \dot{v_2} \\ \dot{a_2} \end{bmatrix} = \begin{bmatrix} v_1\\ a_1 \\ -2v_2^2 + x_2 \\ v_2 \\ a_2 \\ -a_1^3 + v_2 + x_1 + \sin(t) \end{bmatrix}$$

A Final Note: Be Careful

Not all solvers work for all ODEs. For example, the ODEs for 3 body motion

Require the DOP853 solver with low values for rtol and atol to solve

Always review the literature on your specific ODE to find an appropriate solver. In addition, play around with the parameters rtol and atol to ensure you get a proper solution

Initial Conditions