Solving for Möbius transformation coefficients

Möbius transformations are functions of the form

f(z) = \frac{az + b}{cz + d}

where adbc ≠ 0.

A Möbius transformation is uniquely determined by its values at three points. Last year I wrote a post that mentioned how to determine the coefficients of a Möbius transformation. There I said

The unique bilinear transform sending z1, z2, and z3 to w1, w2, and w3 is given by

\frac{(w - w_2)(w_3 - w_1)}{(w - w_1)(w_3 - w_2)} = \frac{(z - z_2)(z_3 - z_1)}{(z - z_1)(z_3 - z_2)}

Plug in your constants and solve for w = f(z).

This is correct, but it still leaves a bit of work to do. In a particular case it’s not hard to find the coefficients, but it would be harder to find the coefficients in general.

There is an explicit formula for each of the parameters a, b, c, and d given the specified points z1, z2, z3 and their images w1, w2, w3 which we will present shortly. We could easily code-up these formulas except for one complication: we may want one of our inputs our outputs to be ∞. This is not merely an edge case: in applications, say to signal processing, you often want to specify the location of poles.

Simplest case: no infinities

If none of the inputs or outputs are infinite, the coefficients are given by

\begin{align*} a &= \begin{vmatrix} z_1w_1 & w_1 & 1 \\ z_2w_2 & w_2 & 1 \\ z_3w_3 & w_3 & 1 \end{vmatrix} \\ b &= \begin{vmatrix} z_1w_1 & z_1 & w_1 \\ z_2w_2 & z_2 & w_2 \\ z_3w_3 & z_3 & w_3 \end{vmatrix} \\ c &= \begin{vmatrix} z_1\phantom{w_1} & w_1 & 1 \\ z_2\phantom{w_1} & w_2 & 1 \\ z_3\phantom{w_1} & w_3 & 1 \end{vmatrix} \\ d &= \begin{vmatrix} z_1w_1 & z_1 & 1 \\ z_2w_2 & z_2 & 1 \\ z_3w_3 & z_3 & 1 \end{vmatrix} \end{align*}

This is one solution; multiplying all four coefficients by a non-zero constant gives another solution. But modulo a constant the solution is unique. To put it another way, the Möbius transformation is unique but it’s representation is only unique up to a constant multiplying the numerator and denominator.

You might hope for a second that this formula would work if you just let floating point infinities take care of themselves. But if any of the z‘s or w‘s is infinite, every determinant above is infinite.

Dealing with infinities

The possible cases of infinities are:

  1. No infinities
  2. Only one z is infinite
  3. Only one w is infinite
  4. A z and a w with the same subscript are infinite
  5. A z and a w with different subscripts are infinite

The order of our z‘s and w‘s is arbitrary, so we can rearrange them if necessary for our convenience. So without loss of generality we may assume

  1. No infinities
  2. Only z1 is infinite
  3. Only w1 is infinite
  4. z1 = w1 = ∞
  5. z1 = w2 = ∞

To handle the case (2) of z1= ∞, divide each of the equations above by z1 and take the limit as z1 approaches ∞. Since dividing one row of a matrix by a constant divides its determinant by the same amount, in each case we divide the first row by z1. This works out nicely because z1 only ever appears in the first row of each matrix. We can handle the case (3) of 21= ∞ analogously.

To handle the case (4) of z1 = w1 = ∞ we divide the first rows by z1 w1 and take the limit. To handle the case (5) of z1 = w2 = ∞ we divide the first rows by z1 and the second rows by w1 and take a limit.

If you find this limiting business dubious, it doesn’t matter: if the result is correct, it’s correct. And since Möbius transformations are determined by their values at three points, you can verify that each case is correct by sticking in z1, z2, z3 and checking that you get w1, w2, w3 out. You could do this manually, which I have, or trust the output of the code at the bottom of the post.

So here are our solutions.

(1) No infinities

Given above.

(2) Only z1 = ∞

a = w1 (w2w3)
b = w1 (z2 w3z3 w2) + w2 w3 (z3z2)
c = w2w3
d = w
1 (z2z3)z2 w2 + z3 w3

(3) Only w1 = ∞

a = z1 (w2 – w3) – z2 w2 + z3 w3
b = z1 (z2 w3z3 w2) + z2 z3 (w2w3)
c = z3 – z2
d = z1 (z2z3)

NB: You could derive these from the case z1 = ∞ by inverting the transform. This means swapping z‘s with w‘s, swapping a and d, and negating b and c.

(4) z1 = w1 = ∞

a = w2w3
b = z2 w3z3 w2
c = 0
d = z2z3

(5) z1 = w2 = ∞

a = w1
b = –z2 w3 + z3 (w3w1)
c = 1
d = –z2

Python code

import numpy as np

def all_finite(z1, z2, z3, w1, w2, w3):
    a = np.linalg.det(
            [[z1*w1, w1, 1],
             [z2*w2, w2, 1],
             [z3*w3, w3, 1]])
    b = np.linalg.det(
            [[z1*w1, z1, w1],
             [z2*w2, z2, w2],
             [z3*w3, z3, w3]])
    c = np.linalg.det(
            [[z1, w1, 1],
             [z2, w2, 1],
             [z3, w3, 1]])
    d = np.linalg.det(
            [[z1*w1, z1, 1],
             [z2*w2, z2, 1],
             [z3*w3, z3, 1]])
    return (a, b, c, d)

def z1_infinite(z1, z2, z3, w1, w2, w3):
    assert(np.isinf(z1))
    a = w1*(w2 - w3)
    b = w1*(z2*w3 - z3*w2) + w2*w3*(z3 - z2)
    c = w2 - w3
    d = w1*(z2 - z3) - z2*w2 + z3*w3
    return (a, b, c, d)    

def w1_infinite(z1, z2, z3, w1, w2, w3):
    assert(np.isinf(w1))
    a = z1*(w2 - w3) - z2*w2 + z3*w3
    b = z1*(z2*w3 - z3*w2) + z2*z3*(w2 - w3)
    c = z3 - z2
    d = z1*(z2 - z3)
    return (a, b, c, d)        

def z1w1_infinite(z1, z2, z3, w1, w2, w3):
    assert(np.isinf(z1) and np.isinf(w1))
    a = w2 - w3
    b = z2*w3 - z3*w2
    c = 0
    d = z2 - z3
    return (a, b, c, d)            

def z1w2_infinite(z1, z2, z3, w1, w2, w3):
    assert(np.isinf(z1) and np.isinf(w2))
    a = w1
    b = -z2*w3 + z3*(w3 - w1)
    c = 1
    d = -z2
    return (a, b, c, d)                

def mobius_coeff(z1, z2, z3, w1, w2, w3):

    infz = np.isinf(z1) or np.isinf(z2) or np.isinf(z3)
    infw = np.isinf(w1) or np.isinf(w2) or np.isinf(w3)    

    if infz:
        if np.isinf(z2):
            z1, z2 = z2, z1
            w1, w2 = w2, w1
        if np.isinf(z3):
            z1, z3 = z3, z1
            w1, w3 = w3, w1
        if infw:
            if np.isinf(w1):
                return z1w1_infinite(z1, z2, z3, w1, w2, w3)
            if np.isinf(w3):
                z2, z3 = z3, z2
                w2, w3 = w3, w2
            return z1w2_infinite(z1, z2, z3, w1, w2, w3)
        else:
            return z1_infinite(z1, z2, z3, w1, w2, w3)
        
    if infw: # and all z finite
        if np.isinf(w2):
            z1, z2 = z2, z1
            w1, w2 = w2, w1
        if np.isinf(w3):
            z1, z3 = z3, z1
            w1, w3 = w3, w1
        return w1_infinite(z1, z2, z3, w1, w2, w3)

    return all_finite(z1, z2, z3, w1, w2, w3)

def mobius(x, a, b, c, d):
    if np.isinf(x):
        if c == 0:
            return np.inf
        return a/c
    if c*x + d == 0:
        return np.inf
    else:
        return (a*x + b)/(c*x + d)
    
def test_mobius(z1, z2, z3, w1, w2, w3):
    tolerance = 1e-6
    a, b, c, d = mobius_coeff(z1, z2, z3, w1, w2, w3)
    for (x, y) in [(z1, w1), (z2, w2), (z3, w3)]:
        m = mobius(x, a, b, c, d)
        assert(np.isinf(m) and np.isinf(y) or abs(m - y) <= tolerance)

test_mobius(1, 2, 3, 6, 4, 2)
test_mobius(1, 2j, 3+7j, 6j, -4, 2) 
test_mobius(np.inf, 2, 3, 8j, -2, 0)
test_mobius(0, np.inf, 2, 3, 8j, -2)
test_mobius(0, -1, np.inf, 2, 8j, -2)
test_mobius(1, 2, 3, np.inf, 44j, 0)
test_mobius(1, 2, 3, 1, np.inf, 40j)
test_mobius(-1, 0, 3j, 1, -1j, np.inf)
test_mobius(np.inf, -1j, 5, np.inf, 2, 8)
test_mobius(1, np.inf, -1j, 5, np.inf, 2)
test_mobius(12, 0, np.inf, -1j, 5, np.inf)
test_mobius(np.inf, -1j, 5, 0, np.inf, -1)
test_mobius(6, np.inf, -1j, 0, 8, np.inf)
test_mobius(6, 3j, np.inf, -1j, np.inf, 1)

More Möbius transformation posts

2 thoughts on “Solving for Möbius transformation coefficients

  1. Using the scipy.linalg.null_space command can make this a slightly more concise.
    Each pair of points corresponds to an equation of the form
    z a + b – z w c – w d = 0.
    We divide through by z and/or w if they are infinite.
    I did have to relax the test, so it does not give exactly the same kind of solution.
    Not sure if posting code works.

    import numpy as np
    import scipy.linalg as sla
    
    # w = (a*z+b)/(c*z+d)
    # z a + b - z w c - w d = 0
    
    def mobius_row(z, w):
        if np.isinf(z):
            if np.isinf(w):
                row = [0, 0, -1, 0]
            else:
                row = [1, 0, -w, 0]
        else:
            if np.isinf(w):
                row = [0, 0, -z, -1]
            else:
                row = [z, 1, -z*w, -w]
        return row
    
    def mobius_coeff(z1, z2, z3, w1, w2, w3):
        A = np.array([ mobius_row(z1, w1),
                       mobius_row(z2, w2),
                       mobius_row(z3, w3) ])
        a, b, c, d = sla.null_space(A).flat
        return a, b, c, d
    
    
    assert(np.isinf(m) and np.isinf(y) or
           np.isinf(m) and 1/abs(y) <= tolerance or
           np.isinf(y) and 1/abs(m) <= tolerance or
           abs(m - y) <= tolerance)
    
  2. Thanks. That’s neat.

    I edited your comment to turn <code> into <pre> and that preserved your code formatting.

Leave a Reply

Your email address will not be published.