Ellipse Construction Advice
I propose an adaptive approach :

The idea is to start from the simple x = a * cos(t) and y = b * sin(t).

Then calculate all the distances between the consecutive points, and progressively adapt the angle between them relatively to each distance compared to the average distance.

Note : this is to improve, as the last point (on top in the image below) does not "move". So its own angle does not vary as the others... edit : this note is false. After some tests, this is just matter of setting the good precision to reach (for instance set the precision to 0.00001 and all is ok and needs 6 iterations, so 6 * the amount of vertices calculations). 

import bpy, bmesh
from math import radians, sqrt, cos, sin

rad90 = radians( 90.0 )
rad180 = radians( 180.0 )

def createVertex( bm, x, y ): #uses bmesh to create a vertex
    return bm.verts.new( [x, y, 0] )

def listSum( list, index ): #helper to sum on a list
    sum = 0
    for i in list:
        sum = sum + i[index]
    return sum

def calcLength( points ): #calculate the lenghts for consecutives points
    prevPoint = points[0]
    for point in points :
        dx = point[0] - prevPoint[0]
        dy = point[1] - prevPoint[1]
        dist = sqrt( dx * dx + dy *dy )
        point[3] = dist
        prevPoint = point

def calcPos( points, a, b ): #calculate the positions following the angles
    angle = 0
    for i in range( 1, len(points) - 1 ):
        point = points[i]
        angle += point[2]
        point[0] = a * cos( angle )
        point[1] = b * sin( angle )

def adjust( points ): #adjust the angle by comparing each length to the mean length
    totalLength = listSum( points, 3 )
    averageLength = totalLength / (len(points) - 1)

    maxRatio = 0
    for i in range( 1, len(points) ):
        point = points[i]
        ratio = (averageLength - point[3]) / averageLength
        point[2] = (1.0 + ratio) * point[2]
        absRatio = abs( ratio )
        if absRatio &gt; maxRatio:
            maxRatio = absRatio
    return maxRatio

def ellipse( bm, a, b, steps, limit ):

    delta = rad90 / steps

    angle = 0.0

    points = [] #will be a list of [ [x, y, angle, length], ...]
    for step in range( steps  + 1 ) :
        x = a * cos( angle )
        y = b * sin( angle )
        points.append( [x, y, delta, 0.0] )
        angle += delta

    doContinue = True
    while doContinue:
        calcLength( points )
        maxRatio = adjust( points )
        calcPos( points, a, b )

        doContinue = maxRatio &gt; limit
        print( maxRatio )

    verts = []    #edit : to add the edges
    for point in points:
        verts.append( createVertex( bm, point[0], point[1] ) )

    for i in range( 1, len(verts) ):
        bm.edges.new( [verts[i - 1], verts[i]] )


A = 2
B = 1

bm = bmesh.new()

ellipse( bm, A, B, 12, 0.01 )

mesh = bpy.context.object.data      
bm.to_mesh(mesh)
mesh.update()


The algorithm seems to converge (stop) around 4 iterations for a 1% around the mean distance and around 6 iterations to be under 0.01%.

Not very mathematical, but it avoids complex analytical calculations...

 (edited : added edges)