Starting to understand Polar Decomposition

PolarDecomp.fs

// Trigonometry
let pi = System.Math.PI
let sin t = System.Math.Sin(t)
let cos t = System.Math.Cos(t)
let asin r = System.Math.Asin(r)
let acos r = System.Math.Acos(r)

// Linear Algebra
type Matrix2x2 = { m11: double; m12: double;
                   m21: double; m22: double; }

let MakeMatrix2x2 a b c d = { m11=a; m12=b;
                              m21=c; m22=d; }
let MakeRotMatrix t = MakeMatrix2x2 (cos t) (-sin t) (sin t) (cos t)
let MakeScaleMatrix x y = MakeMatrix2x2 x 0.0 0.0 y

let Mult2x2 a b = { m11 = a.m11*b.m11 + a.m12*b.m21;
                    m12 = a.m11*b.m12 + a.m12*b.m22;
                    m21 = a.m21*b.m11 + a.m22*b.m21;
                    m22 = a.m21*b.m12 + a.m22*b.m22; }
let Add2x2 a b = { m11=a.m11+b.m11; m12=a.m12+b.m12;
                   m21=a.m21+b.m21; m22=a.m22+b.m22; }
let Scale2x2 s a = { m11=s*a.m11; m12=s*a.m12;
                     m21=s*a.m21; m22=s*a.m22; }
let Transpose2x2 a = { a with m12=a.m21; m21=a.m12; }
let Det2x2 a = a.m11*a.m22-a.m12*a.m21
let Inverse2x2 a = Scale2x2 (1.0/Det2x2 a) { m11=    a.m22; m12=0.0-a.m12;
                                             m21=0.0-a.m21; m22=    a.m11; }
let Print2x2 a = printfn "[ [ %.6f, %.6f ], [ %.6f, %.6f ] ]" a.m11 a.m12 a.m21 a.m22

let Decomp2x2 a =
  // Iteratively average the matrix entrywise with the transpose of its inverse. (!)
  // Thanks: (Shoemake 92), (Highham 86), etc.
  let q = ref a
  for i=1 to 5 do
    q := Scale2x2 0.5 (Add2x2 !q (Transpose2x2 (Inverse2x2 !q)))
  (!q, Mult2x2 (Inverse2x2 !q) a)

// Test Polar Decomp
let t = 2.0*pi/17.0
let sx = 2.0
let sy = 3.0
printfn "Actual: t=%.6f, sx=%.6f, sy=%.6f" t sx sy

let r = MakeRotMatrix t
let s = MakeScaleMatrix sx sy
let m = Mult2x2 r s

let (rGuess, sGuess) = Decomp2x2 m

let tGuess = 0.25*((acos (rGuess).m11)+(asin -(rGuess).m12)+
                   (asin (rGuess).m21)+(acos  (rGuess).m22))
let sxGuess = sGuess.m11
let syGuess = sGuess.m22
printfn "Decomp: t=%.6f, sx=%.6f, sy=%.6f" tGuess sxGuess syGuess

// Actual: t=0.369599, sx=2.000000, sy=3.000000
// Decomp: t=0.369599, sx=2.000000, sy=3.000000