def Dirichlet_plane(Nx, Ny, x1, y1, x2, y2):
    hx = 1/ Nx
    hy = 1/ Ny

    xnodes = np.linspace(x1, x2, Nx + 1)
    xnodes = np.repeat(xnodes, Ny + 1)
    xnodes = np.reshape(xnodes, (Nx + 1, Ny + 1))

    ynodes = np.linspace(y1, y2, Ny + 1)
    ynodes = np.repeat(ynodes, Nx + 1)
    ynodes = np.reshape(ynodes, (Ny + 1, Nx + 1)).T

    # Заполняем таблицу узлов
    nodes = np.stack((xnodes, ynodes), axis=-1)

    u = np.zeros((Nx + 1, Ny + 1))
    f = np.zeros((Nx + 1, Ny + 1))

    # Заполняем матрицу для функции f
    f[:, :] = fxy(nodes[..., 0], nodes[..., 1])
    f[1, 1:Ny] += gxy(nodes[0, 1:Ny][:, 0], nodes[0, 1:Ny][:, 1]) / hx / hx

    # Заполняем граничные условия для u
    u[0, :] = gxy(nodes[0, :, 0], nodes[0, :, 1])
    u[Nx, :] = gxy(nodes[Nx, :, 0], nodes[Nx, :, 1])
    u[:, 0] = gxy(nodes[:, 0, 0], nodes[:, 0, 1])
    u[:, Ny] = gxy(nodes[:, Ny, 0], nodes[:, Ny, 1])
    print(u)

    # Заполняем матрицу для \phi_k(n)
    phi = ifft(f[1:Nx, 1:Ny].T).T

    # Заполняем матрицу для c_k(n)
    c = np.zeros((Nx - 1, Ny + 1))
    c[:, 0] = ifft(f[1:Nx, 0])
    c[:, Ny] = ifft(f[1:Nx, Ny])

    # Решаем систему уравнений
    for k in range(1, Nx):
        b = np.array(phi[k - 1])
        b[0] += c[k - 1, 0] / hy / hy
        b[Ny - 2] += c[k - 1, Ny] / hy / hy

        A = np.zeros((Ny - 1, Ny - 1))
        flat = np.ravel(A)
        flat[Ny-1::Ny] = -1 / hy / hy
        flat[1::Ny] = -1 / hy / hy
        flat[0::Ny] = 2 / hy / hy + 4 * (np.sin(math.pi * k * hx / 2) ** 2)
        c[k - 1, 1:Ny] = np.linalg.solve(A, b)

    # Находим функцию u
    u[1:Nx, 1:Ny] = ifft(c[:, 1:Ny].T).T
    return u