def _get_elem_neighbors(p, nelx, nely, nelz):
    """
    Gets the elements which surround the element in point p1 (no diagonals).
    """

    eset = {}
    if nelz == 0:
        e1 = (max(p[0]-1, 0), p[1])
        e2 = (min(p[0]+1, nely-1), p[1])
        e3 = (p[0], max(p[1]-1, 0))
        e4 = (p[0], min(p[1]+1, nelx-1))

        eset = {e1, e2, e3, e4}
    else:
        e1 = (max(p[0]-1, 0), p[1], p[2])
        e2 = (min(p[0]+1, nelz-1), p[1], p[2])
        e3 = (p[0], max(p[1]-1, 0), p[2])
        e4 = (p[0], min(p[1]+1, nely-1), p[2])
        e5 = (p[0], p[1], max(p[2]-1, 0))
        e6 = (p[0], p[1], min(p[2]+1, nelx-1))

        eset = {e1, e2, e3, e4, e5, e6}

    return list(eset)