Changeset 3820


Ignore:
Timestamp:
Feb 16, 2015, 10:21:01 AM (7 years ago)
Author:
sam
Message:

math: disable unstable LU decomposition for matrix inversion.

Location:
trunk/src
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/lol/math/matrix.h

    r3807 r3820  
    548548T determinant(mat_t<T, N, N> const &m)
    549549{
     550#if 0 /* XXX: LU decomposition is currently broken */
    550551    mat_t<T, N, N> L, U;
    551552    lu_decomposition(m, L, U);
     
    556557
    557558    return det;
     559#else
     560    T ret = (T)0;
     561    for (int i = 0; i < N; ++i)
     562        ret += m[i][0] * cofactor(m, i, 0);
     563    return ret;
     564#endif
     565}
     566
     567template<typename T>
     568T const & determinant(mat_t<T, 1, 1> const &m)
     569{
     570    return m[0][0];
    558571}
    559572
     
    615628mat_t<T, N, N> inverse(mat_t<T, N, N> const &m)
    616629{
     630#if 0 /* XXX: LU decomposition is currently broken */
    617631    mat_t<T, N, N> L, U;
    618632    lu_decomposition(m, L, U);
     
    621635
    622636    return ret;
     637#else
     638    mat_t<T, N, N> ret;
     639
     640    T d = determinant(m);
     641    if (d)
     642    {
     643        /* Transpose result matrix on the fly */
     644        for (int i = 0; i < N; ++i)
     645            for (int j = 0; j < N; ++j)
     646                ret[i][j] = cofactor(m, j, i) / d;
     647    }
     648    return ret;
     649#endif
    623650}
    624651
  • trunk/src/t/math/matrix.cpp

    r3817 r3820  
    7777    lolunit_declare_test(inverse_2x2)
    7878    {
    79         mat2 m0 = inv2;
    80         mat2 m1 = inverse(m0);
    81         mat2 m2 = m0 * m1;
    82 
     79        mat2 m(vec2(4, 3),
     80               vec2(3, 2));
     81
     82        mat2 m1 = inverse(m);
     83        for (int j = 0; j < 2; ++j)
     84        for (int i = 0; i < 2; ++i)
     85            lolunit_assert_equal(m1[i][j], m1[i][j]);
     86
     87        mat2 m2 = m1 * m;
    8388        for (int j = 0; j < 2; ++j)
    8489        for (int i = 0; i < 2; ++i)
     
    8691    }
    8792
     93#if 0 /* XXX: LU decomposition is currently broken */
    8894    lolunit_declare_test(lu_decomposition_3x3)
    8995    {
     
    227233            lolunit_assert_doubles_equal(m2[i][j], mat4(1.f)[i][j], 1e-5);
    228234    }
     235#endif
    229236
    230237    lolunit_declare_test(inverse_3x3)
     
    233240               vec3(3, 2, 3),
    234241               vec3(9, 5, 7));
    235         mat3 m2 = inverse(m) * m;
    236 
     242
     243        mat3 m1 = inverse(m);
     244        for (int j = 0; j < 3; ++j)
     245        for (int i = 0; i < 3; ++i)
     246            lolunit_assert_equal(m1[i][j], m1[i][j]);
     247
     248        mat3 m2 = m1 * m;
    237249        for (int j = 0; j < 3; ++j)
    238250        for (int i = 0; i < 3; ++i)
     
    246258               vec4( 4,  2,  5, -4),
    247259               vec4( 5, -3, -7, -6));
    248         mat4 m2 = inverse(m) * m;
    249 
     260
     261        mat4 m1 = inverse(m);
     262        for (int j = 0; j < 4; ++j)
     263        for (int i = 0; i < 4; ++i)
     264            lolunit_assert_equal(m1[i][j], m1[i][j]);
     265
     266        mat4 m2 = m1 * m;
    250267        for (int j = 0; j < 4; ++j)
    251268        for (int i = 0; i < 4; ++i)
     
    259276               vec4(0, -1,  0,  0),
    260277               vec4(0,  0, -1,  1));
    261         mat4 m2 = inverse(m) * m;
    262 
     278        mat4 m1 = inverse(m);
     279        for (int j = 0; j < 4; ++j)
     280        for (int i = 0; i < 4; ++i)
     281            lolunit_assert_equal(m1[i][j], m1[i][j]);
     282
     283        mat4 m2 = m1 * m;
    263284        for (int j = 0; j < 4; ++j)
    264285        for (int i = 0; i < 4; ++i)
Note: See TracChangeset for help on using the changeset viewer.