Kinematic model for a 5axis mill with universal (nutating) head

I would first try using the exact algebraic expression (ie sin(50°), cos(50°)) and see what the rotation matrix looks like.

If you manage to get the correct forward and inverse kinematic model using pen and paper then my hat is definitely going to be off to you. I started out that way ...
I attach my jupyter notebook (remove the '.txt') which contains all the ipython symbolic math to automatically derive all the formulas in the documentation. Using jupyter and ipython will save you a LOT of paper and headache.
You will need to adjust some of the matrices and order of the matrix multiplications though.

I didn't bother to type out that part of the documentation and my hand notes are long gone.
I'll see if I can retrace my steps but it may take me a few days to do so.

Ok, now it all came back to me:

To calculate the joint angles (which is done in the function 'kins_calc_jnt_angles' in 'remap.py' ) we use the rotation matrix part (ie disregard the fourth row with all the zeros ) of the transformation matrix as derived for the inverse tool kinematic:

 

the columns of this matrix represent the tool orientation vectors [tx, ty, tz] and hence the 3rd column represents the tool-z orientation vector we are trying to orient the spindle head to.
So using the variable names from the remap we get:

Kwx = Cb*t + Sb * u
Kwy = v
Kwz = -Sb*t + Cb*u

with the substitution for 'v' we can then solve for Sa and thus for the angle of joint A 
Then we use the two remaining equations to solve for Sb and thus for the angle of Joint B

The output of this calculation needs to conditioned so the correct angle values are output for the machine rotary joints.

As a side note:

in the function ' kins_calc_tool_rot_c_for_horizontal_x' of remap.py there is  this comment:

    # Mathematically we take the symbolic formula found in row 3, column 1 of the transformation
    # matrix from the inverse tool-kinematics, equal that to zero and solve for 'tc'.
    # ie: -(tctc-vstc)cb-(cactc-tstc)sb=0

  Note that the formula given in the comment may not be readily recognizable in the matrix shown above yet it is indeed identical to the third row of the first column. It seems that in the final version of the documentation 'Ctc' and 'Stc' were factored out instead of 'Cb' and 'Sb'.

 

Correction:
with the substitution for 'v' we can then solve for Sa and thus for the angle of joint A 

with the substitution for 'v' we can then solve for Ca and thus for the angle of joint A