Re-arranging verticies
Displace a grid.

By nature  moving an imported obj by the noise result is going to look voxelated.
Could instead simply displace the matching z coordinates of a grid by the noise. (Akin to using as a displacement texture)
import numpy as np
import bpy

def interpolant(t):
    return t * t * t * (t * (t * 6 - 15) + 10)

def generate_perlin_noise_2d(
        shape, res, tileable=(False, False), interpolant=interpolant
):

    delta = (res[0] / shape[0], res[1] / shape[1])
    d = (shape[0] // res[0], shape[1] // res[1])
    grid = np.mgrid[0:res[0]:delta[0], 0:res[1]:delta[1]] \
               .transpose(1, 2, 0) % 1
    # Gradients
    angles = 2 * np.pi * np.random.rand(res[0] + 1, res[1] + 1)
    gradients = np.dstack((np.cos(angles), np.sin(angles)))
    if tileable[0]:
        gradients[-1, :] = gradients[0, :]
    if tileable[1]:
        gradients[:, -1] = gradients[:, 0]
    gradients = gradients.repeat(d[0], 0).repeat(d[1], 1)
    g00 = gradients[:-d[0], :-d[1]]
    g10 = gradients[d[0]:, :-d[1]]
    g01 = gradients[:-d[0], d[1]:]
    g11 = gradients[d[0]:, d[1]:]
    # Ramps
    n00 = np.sum(np.dstack((grid[:, :, 0], grid[:, :, 1])) * g00, 2)
    n10 = np.sum(np.dstack((grid[:, :, 0] - 1, grid[:, :, 1])) * g10, 2)
    n01 = np.sum(np.dstack((grid[:, :, 0], grid[:, :, 1] - 1)) * g01, 2)
    n11 = np.sum(np.dstack((grid[:, :, 0] - 1, grid[:, :, 1] - 1)) * g11, 2)
    # Interpolation
    t = interpolant(grid)
    n0 = n00 * (1 - t[:, :, 0]) + t[:, :, 0] * n10
    n1 = n01 * (1 - t[:, :, 0]) + t[:, :, 0] * n11
    return np.sqrt(2) * ((1 - t[:, :, 1]) * n0 + t[:, :, 1] * n1)

def generate_fractal_noise_2d(
        shape, res, octaves=1, persistence=0.5,
        lacunarity=2, tileable=(False, False),
        interpolant=interpolant
):
    noise = np.zeros(shape)
    frequency = 1
    amplitude = 1
    for _ in range(octaves):
        noise += amplitude * generate_perlin_noise_2d(
            shape, (frequency * res[0], frequency * res[1]), tileable, interpolant
        )
        frequency *= lacunarity
        amplitude *= persistence
    return noise
# lightweight test
#p = generate_fractal_noise_2d((16, 16), (2,  2), 3)
p = generate_fractal_noise_2d((256, 256), (8, 8), 5)
print(&quot;-&quot; * 33)
i, j = p.shape
scale = 0.5

bpy.ops.mesh.primitive_grid_add(
        x_subdivisions=i,
        y_subdivisions=j,
        size=10,
        )
        
grid = bpy.context.object

for v in grid.data.vertices:
    v.co.z = scale * p[v.index // i][v.index % j]

Please note, have used a loop by way of example, in a rush (off to the vet)  and will get back
Far quicker to  use
grid.data.vertices.foreach_set(&quot;co&quot;, flat_array)

and set with  a ravelled list of coordinates
[x0, y0, z0, x1, y1, z1, ....]

eg making x and y with linspace or meshgrid, then slicing in z as shown in https://blender.stackexchange.com/a/172973/15543 is according to stackoverflow link quickest way to ravel.
Note blender has a perlin noise generator as well, will also try and find link re this.