#216 Inverse of a matrix and matrix determinant routines

[Nyah Check]

This is a new implementation of the bn_mat_inverse() and bn_mat_determinant routines in src/libbn/mat.c of BRL-CAD; with the following amendments:

bn_mat_inverse:

• The matrix pointed at by "input" is inverted and stored in the area
• pointed at by "output".
  *
• Invert a 4-by-4 matrix using direct computation.
• Uses 100 multiplications and 61 additions or substractions, total

bn_mat_determinant:

• Calculates the determinant of the 4X4 matrix
• (This implementation requires only 34 multiplications
• and 20 additions or substractions, saving 6 multiplications
• and 3 additions or substractions compared to the previous one.

[Nyah Check]

Here are the test results comparing the old bn_mat_determinant to the new:

Old_determinant Test:

Det: -3291.00
real 1m5.151s
user 1m5.081s
sys 0m0.019s

New_determinant Test:

Det: -3291.00
real 0m54.921s
user 0m54.815s
sys 0m0.054s

[Nyah Check]

here is another summary:
Summary:

bn_mat_determinant()

GCC 4.6.3 compiles both existing libbn and my implementations to roughly equally efficient binary on x86-64 if optimizations are enabled. Without compiler optimizations, mine is clearly faster.

Without optimizations, mine runs in 154 cycles (minimum and median), whereas the libbn one takes 210 cycles minimum, 219 cycles median. (Mine is almost a third faster!)

With typical optimizations, '-O2', mine runs in 92 or 93 cycles (minimum and median), whereas the libbn one runs in 95 cycles (minimum and median).

The relative difference in the results for a large set of random matrices was always under one billionth, 10-9. That is, the error is roughly the same order as the precision of a square root of some value, sqrt().

bn_mat_inverse()

Using GCC-4.6.3, my implementation is demonstrably faster than the existing libbn implementation -- even if the existing libbn implementation cannot invert all non-singular real-valued 4x4 matrices! (I verified some of my random matrices using Maple.)

Without optimization, my implementation is three times as fast, completing in 544 cycles (median value), whereas libbn implementation requires 1660 cycles (median value).

With -O2, my implementation takes 423 (minimum) to 432 (median) cycles, whereas libbn implementation takes 610 (minimum) to 660 (median) cycles.

My preferred maximum level of optimization is '-O3 -fomit-frame-pointer'; my implementation takes 423 (minimum) to 427 (median) cycles, whereas libbn implementation takes 466 (minimum) to 483 (median) cycles; in the median case, my implementation is still 11% faster.

The maximum level of optimization I tried is '-O3 -fomit-frame-pointer -ftree-vectorize -march=native -mtune=native'. Mine still runs in 420 to 426 cycles, whereas libbn took 480 to 500 cycles.

[Nyah Check]

Wrote Tests to compare Inverse of Old matrix inverse routine and new_matrix Inverse routine which returned the following results:

Old Inverse matrix Test:

real 8m7.388s
user 8m6.620s
sys 0m0.223s

New Inverse matrix Test:

real 4m13.705s
user 4m12.766s
sys 0m0.212s

[Nyah Check]