Revision history [back]

The following alternative code is running without errors. The Möbius transform is slightly changed, and the $\sqrt 2=|s|$ was added. (I think it is also wanted.)

K.<s> = QuadraticField( -2 )
OK    = K.ring_of_integers()
IP    = ProjectiveSpace( 1, K )

points = ( -1, 0, 1 )
def lift_point( x ):
    return vector( [ x, 1 ] )

def FLT( M, z ):
    """takes the fractional linear transformation/Mobius transformation of a z in CC"""
    Mz = IP( ( M * lift_point(z) ).list() )
    if Mz == IP( (1,0) ):
        return 'infinity'
    return Mz[0]

def components( z ):
    return ( z[0], z[1]*sqrt(2) )

var('a','b','r','x','y');

func_cir = (x-a)**2 + (y-b)**2 == r**2

def get_circles( number_of_circles=10, show=False ):
    j = 0
    circles  = []
    while j < number_of_circles:

        M = random_matrix( OK, 2, algorithm = 'unimodular' )
        A, B, C = [ components( FLT( M, z ) ) for z in points ]

        print ( "\n%s\nM is %s\nA,B,C are %s, %s, %s ."
                % ( 30*'=', M.list(), A, B, C ) )
        colinearityMatrix = matrix( 3,3, [ A[0], A[1], 1 ,
                                           B[0], B[1], 1 ,
                                           C[0], C[1], 1 ] )
        colinearity = bool( colinearityMatrix.det() == 0 )
        # print ( "Are A,B,C colinear? %s" % colinearity )
        if colinearity:
            print "!!! COLINEAR... IGNORED !!!\n"

        elif 'infinity' in [ A, B, C ]:
            print "!!! INFINITY... IGNORED !!!\n"

        else:
            j += 1

            eq1 = func_cir.subs( x==A[0] ).subs( y==A[1] )
            eq2 = func_cir.subs( x==B[0] ).subs( y==B[1] )
            eq3 = func_cir.subs( x==C[0] ).subs( y==C[1] )

            sol = solve( [ eq1, eq2, eq3 ], a,b,r )

            print "%s\n[%s]\nSolutions:\n%s" % ( 30*'-', j, sol )
            c = circle( ( sol[0][0].rhs(),
                          sol[0][1].rhs() ),
                          abs( sol[0][2].rhs() ) )
            circles.append( c )

    if show:
        plot_object = plot( [] )
        for c in circles:
            plot_object += plot( c )
        plot_object.show()
    return circles

Now the following calls did the job for me without errors:

get_circles( 1000 );

There is a verbose, aggressive print show in between, but the message is, there are no errors.

Let us also generate three circles, and show them:

get_circles( 3, show=True  );

The following alternative code is running without errors. The Möbius transform is slightly changed, and the $\sqrt 2=|s|$ was added. (I think it is also wanted.)

K.<s> = QuadraticField( -2 )
OK    = K.ring_of_integers()
IP    = ProjectiveSpace( 1, K )

points = ( -1, 0, 1 )
def lift_point( x ):
    return vector( [ x, 1 ] )

def FLT( M, z ):
    """takes the fractional linear transformation/Mobius transformation of a z in CC"""
    Mz = IP( ( M * lift_point(z) ).list() )
    if Mz == IP( (1,0) ):
        return 'infinity'
    return Mz[0]

def components( z ):
    if z == 'infinity':
        return None
    return ( z[0], z[1]*sqrt(2) )

var('a','b','r','x','y');

func_cir = (x-a)**2 + (y-b)**2 == r**2

def get_circles( number_of_circles=10, show=False ):
    j = 0
    circles  = []
    while j < number_of_circles:

        M = random_matrix( OK, 2, algorithm = 'unimodular' )
        A, B, C = [ components( FLT( M, z ) ) for z in points ]

        print ( "\n%s\nM is %s\nA,B,C are %s, %s, %s ."
                % ( 30*'=', M.list(), A, B, C ) )
        noneOccured = bool( None in ( A, B, C ) ) 
        colinearityMatrix = matrix( 3,3, [ A[0], A[1], 1 ,
                                           B[0], B[1], 1 ,
                                           C[0], C[1], 1 ] )
        colinearity = ) if not noneOccured else None

        if noneOccured:
            print "!!! INFINITY... IGNORED !!!\n"

        elif bool( colinearityMatrix.det() == 0 )
):
            # print ( "Are A,B,C colinear? %s" % colinearity )
        if colinearity:
            print "!!! COLINEAR... IGNORED !!!\n"

        elif 'infinity' in [ A, B, C ]:
            print "!!! INFINITY... IGNORED !!!\n"

        else:
            j += 1

            eq1 = func_cir.subs( x==A[0] ).subs( y==A[1] )
            eq2 = func_cir.subs( x==B[0] ).subs( y==B[1] )
            eq3 = func_cir.subs( x==C[0] ).subs( y==C[1] )

            sol = solve( [ eq1, eq2, eq3 ], a,b,r )

            print "%s\n[%s]\nSolutions:\n%s" % ( 30*'-', j, sol )
            c = circle( ( sol[0][0].rhs(),
                          sol[0][1].rhs() ),
                          abs( sol[0][2].rhs() ) )
            circles.append( c )

    if show:
        plot_object = plot( [] )
        for c in circles:
            plot_object += plot( c )
        plot_object.show()
    return circles

Now the following calls did the job for me without errors:

get_circles( 1000 );

There is a verbose, aggressive print show in between, but the message is, there are no errors.

Let us also generate three circles, and show them:

get_circles( 3, show=True  );