Function to compute the discrete 2D Laplacian operator in the frequency domain.

Laplacian operator in frequency domain

import math
import numpy as np
def kv2(Nx, Ny, Lx, Ly):
    """ Compute the discrete 2D Laplacian in the frequency domain."""
    kv2 = np.empty((Nx, Ny))
    # We don't want floor division with Lx, Ly
    Lx = float(Lx)
    Ly = float(Ly)

    def set_point(kx, ky):
        return -4 * math.pi**2 * ((kx/Lx)**2 + (ky/Ly)**2)

    for kx, ky in np.ndindex(Nx, Ny):
        if (kx <= Nx//2) and (ky <= Ny//2):
            kv2[kx, ky] = set_point(kx, ky)
        elif (kx <= Nx//2) and (ky > Ny//2):
            kv2[kx, ky] = set_point(kx, ky - Ny)
        elif (kx > Nx//2) and (ky <= Ny//2):
            kv2[kx, ky] = set_point(kx - Nx, ky)
        else:
            kv2[kx, ky] = set_point(kx - Nx, ky - Ny)

    return kv2