How do I get a cylinder arrangement of arrays with gradient (increasing/decreasing) heights?
You have several flaws in your design.


The normal distribution is usually a two-dimensional function which only takes one variable as an input. You could multiply f(x) and f(y) and use the resulting value as the height.
Nearing a lot of iterations, you shouldn't use bpy.ops
Use the bmesh module to edit meshes and link the data blocks manually to the scene.
Vector objects can be instantiated with a tuple
Blender is slow with a lot of objects, aim to include all geometry in one object in the future




The resulting code may seem convoluted, but is faster and easier to maintain.

import bpy
import bmesh
import math
from mathutils import Vector

RADIUS = 0.16
HEIGHT = 32
SEGMENTS = 32
DISTANCE = 1
COLS = 20
ROWS = 20
SCALE = 1

MEAN = 0
STD_DEVIATION = 4

def normal_dist(x):
    global MEAN, STD_DEVIATION
    mean = MEAN
    std_deviation = STD_DEVIATION

    exponent = - math.pow(x - mean, 2) / (2 * math.pow(std_deviation, 2))
    return 1 / (std_deviation * math.pow(2*math.pi, 0.5)) * math.exp(exponent)

bm = bmesh.new()
bmesh.ops.create_circle(bm, segments=SEGMENTS, diameter=RADIUS*2, cap_ends=True)
bm.faces.ensure_lookup_table()
ret = bmesh.ops.extrude_discrete_faces(bm, faces=bm.faces)
bmesh.ops.translate(bm, vec=Vector((0,0, HEIGHT)), verts=ret['faces'][0].verts)


scn = bpy.context.scene

name = "me"

for i in range(COLS):
    for j in range(ROWS):
        pos = Vector((i - COLS/2, j - COLS/2, 0))
        pos *= DISTANCE

        name = "cylinder_" + str(i) + "_" + str(j)
        if name not in bpy.data.meshes:
            me = bpy.data.meshes.new(name)
        else:
            me = bpy.data.meshes[name]
        bm.to_mesh(me)

        if name not in bpy.data.objects:
            ob = bpy.data.objects.new(name, me)
        else:
            ob = bpy.data.objects[name]
            ob.data = me

        ob.location = pos.copy()
        pos *= SCALE
        ob.scale.z = normal_dist(pos.x) * normal_dist(pos.y) * HEIGHT

        if ob.name not in scn.objects:
            scn.objects.link(ob)


bm.free()



Sooo cool.