Study guide: Finite difference methods for wave motion

$$ \newcommand{\uex}{{u_{\small\mbox{e}}}} \newcommand{\half}{\frac{1}{2}} \newcommand{\halfi}{{1/2}} \newcommand{\xpoint}{\boldsymbol{x}} \newcommand{\normalvec}{\boldsymbol{n}} \newcommand{\Oof}[1]{\mathcal{O}(#1)} \newcommand{\Ix}{\mathcal{I}_x} \newcommand{\Iy}{\mathcal{I}_y} \newcommand{\It}{\mathcal{I}_t} \newcommand{\setb}[1]{#1^0} % set begin \newcommand{\sete}[1]{#1^{-1}} % set end \newcommand{\setl}[1]{#1^-} \newcommand{\setr}[1]{#1^+} \newcommand{\seti}[1]{#1^i} \newcommand{\Real}{\mathbb{R}} $$

Declaring variables and annotating the code

Pure Python code:

def advance_scalar(u, u_1, u_2, f, x, y, t,
                   n, Cx2, Cy2, dt2, D1=2, D2=1):
    Ix = range(0, u.shape[0]);  Iy = range(0, u.shape[1])
    for i in Ix[1:-1]:
        for j in Iy[1:-1]:
            u_xx = u_1[i-1,j] - 2*u_1[i,j] + u_1[i+1,j]
            u_yy = u_1[i,j-1] - 2*u_1[i,j] + u_1[i,j+1]
            u[i,j] = D1*u_1[i,j] - D2*u_2[i,j] + \ 
                     Cx2*u_xx + Cy2*u_yy + dt2*f(x[i], y[j], t[n])

Copy this function and put it in a file with .pyx extension.

Add type of variables:

function(a, b) $ \rightarrow $ cpdef function(int a, double b)

v = 1.2 $ \rightarrow $ cdef double v = 1.2

Array declaration: np.ndarray[np.float64_t, ndim=2, mode='c'] u