def div(p1, p2):
    """Divergence of P=(p1, p2) using backward differece scheme:

    div(P) = p1_x + p2_y

    Parameters
    ----------
    p1 : ~numpy.ndarray
        The first component to be processed.
    p2 : ~numpy.ndarray
        The second component to be processed.

    Returns
    -------
    div_p : ~numpy.ndarray
        The divergence of P=(p1, p2).

    """
    p1_x = p1.copy()
    p1_x[:, :-1] -= p1[:, 1:]
    p2_y = p2.copy()
    p2_y[:-1, :] -= p2[1:, :]
    div_p = p1_x + p2_y
    return div_p
