How do you create a sphere with a certain amount of edges?
Attempt at the D9



The first part of the code, is an algorithm to estimate N equally spaced points on a sphere.

Adapted from 

Constructing Polyhedra from Repelling Points on a Sphere


  by Simon Tatham
  
  Scatters N points randomly on the surface of a sphere, allows them to
  move under mutual inverse-square repulsion while keeping them
  constrained to the sphere's surface, waits until they converge to
  fixed locations, and outputs their locations.


Once the points are determined to within a tolerance TOL they are
rotated to make one the north pole. A 1 radius UV sphere is createdand holes are bisected into the mesh by cutting with a plane normal to
the radial of each point, and moved rat along the 1 unit radius.


Image showing convex hulls of 9 points on left, and sphere with cutting planes on right.

The Code

Copy the code below, paste in text editor and run script. Change the parameters to suit. Note it might hang every now and then finding a solution , if so hit CtrlC in the system console.

import bpy
import bmesh
from mathutils import Vector
import random
from math import pi, asin, atan2, cos, sin, radians

# parmaaters
n = 9  # number of points on sphere
rat = .86 # how far along the radius to bisect
u_segments = 32  # UV sphere settings
v_segments = 32 
thickness = 0.2  # solidify thickness
TOL = 1e-7

points = [Vector((0, 0, 1))]
for i in range(n - 1):
    theta = random.random() * radians(360)
    phi = 2 * asin(random.random() * 2 - 1)
    points.append(Vector((cos(theta) * cos(phi), 
           sin(theta) * cos(phi), 
           sin(phi))))

while True:
    # Determine the total force acting on each point.
    forces = []
    for i in range(len(points)):
        p = points[i]
        f = Vector()
        ftotal = 0
        for j in range(len(points)):
            if j == i: continue
            q = points[j]
            # Find the distance vector, and its length.
            dv = p - q
            dl = dv.length
            dl3 = dl * dl * dl
            fv = dv / dl3
            # Add to the total force on the point p.
            f = f + fv
        # Stick this in the forces array.
        forces.append(f)
        # Add to the running sum of the total forces/distances.
        ftotal = ftotal + f.length

    fscale = 1 if ftotal &lt;= 0.25 else 0.25 / ftotal

    # Move each point, and normalise. While we do this, also track
    # the distance each point ends up moving.
    dist = 0
    for i in range(len(points)):
        p = points[i]
        f = forces[i]
        p2 = (p + fscale * f).normalized()

        dv = p - p2
        dist = dist + dv.length
        points[i] = p2
    # Done. Check for convergence and finish.
    if dist &lt; TOL: # TOL
        break

context = bpy.context
scene = context.scene
# make one point north pole.
R = points[0].rotation_difference(Vector((0, 0, 1))).to_matrix()
points = [R @ p for p in points]

bm = bmesh.new()
#bmesh.ops.create_icosphere(bm, diameter=1, subdivisions=5 )
bmesh.ops.create_uvsphere(bm, 
        diameter=1, 
        u_segments=u_segments, 
        v_segments=v_segments)

for p in points:
    ret = bmesh.ops.bisect_plane(bm, 
            geom=bm.faces[:]+bm.edges[:]+bm.verts[:],
            plane_co= rat * p,
            plane_no=-p,
            clear_outer=False,
            clear_inner=True)
    '''
    # add visible cutting planes (using bisect)
    R = Vector(p).to_track_quat('Z', 'Y').to_matrix().to_4x4()
    R.translation = rat * p
    bmesh.ops.create_grid(bm, x_segments=1, y_segments=1, size=1, matrix=R)
    '''    
me = bpy.data.meshes.new("dice")
bm.to_mesh(me)
ob = bpy.data.objects.new("dice", me)   
scene.collection.objects.link(ob)
context.view_layer.objects.active = ob
ob.select_set(True)
ob.location = scene.cursor.location
# add a thickness modifier
sm = ob.modifiers.new("Solidify", 'SOLIDIFY')
sm.thickness = thickness


Alternatively could hard code in a set of  9 points.

from mathutils import Vector
raw = [
0.798313857169, -0.442003436215, 0.409057389405,
-0.552293193116, -0.798759837392, 0.238652364346,
-0.173114881688, 0.873088163724, -0.455794137855,
0.749969728716, 0.23229964104, -0.619340199553,
-0.872965308015, 0.388085187337, 0.295502044616,
-0.600669921728, -0.124465371314, -0.789749337749,
0.278179133859, -0.797203763248, -0.535800829925,
-0.105061227491, -0.0758822299659, 0.991566450448,
0.477641726235, 0.744841584279, 0.465907067653
]
it = iter(raw)
points = list(map(Vector, zip(it, it, it)))


EDIT:

Updated for 2.8, for prior version  see previous edit revision.