qvm/detail/determinant_impl.hpp

0001 #ifndef BOOST_QVM_DETAIL_DETERMINANT_IMPL_HPP_INCLUDED
0002 #define BOOST_QVM_DETAIL_DETERMINANT_IMPL_HPP_INCLUDED
0003
0004 // Copyright 2008-2022 Emil Dotchevski and Reverge Studios, Inc.
0005
0006 // Distributed under the Boost Software License, Version 1.0. (See accompanying
0007 // file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
0008
0009 #include <boost/qvm/config.hpp>
0010 #include <boost/qvm/mat_traits_array.hpp>
0011 #include <boost/qvm/static_assert.hpp>
0012
0013 namespace boost { namespace qvm {
0014
0015 namespace
0016 qvm_detail
0017     {
0018     template <int N>
0019     struct
0020     det_size
0021         {
0022         };
0023
0024     template <class M>
0025     BOOST_QVM_CONSTEXPR BOOST_QVM_INLINE_TRIVIAL
0026     typename mat_traits<M>::scalar_type
0027     determinant_impl_( M const & a, det_size<2> )
0028         {
0029         return
0030             mat_traits<M>::template read_element<0,0>(a) * mat_traits<M>::template read_element<1,1>(a) -
0031             mat_traits<M>::template read_element<1,0>(a) * mat_traits<M>::template read_element<0,1>(a);
0032         }
0033
0034     template <class M,int N>
0035     BOOST_QVM_CONSTEXPR BOOST_QVM_INLINE_RECURSION
0036     typename mat_traits<M>::scalar_type
0037     determinant_impl_( M const & a, det_size<N> )
0038         {
0039         typedef typename mat_traits<M>::scalar_type T;
0040         T m[N-1][N-1];
0041         T det=T(0);
0042         for( int j1=0; j1!=N; ++j1 )
0043             {
0044             for( int i=1; i!=N; ++i )
0045                 {
0046                 int j2 = 0;
0047                 for( int j=0; j!=N; ++j )
0048                     {
0049                     if( j==j1 )
0050                         continue;
0051                     m[i-1][j2] = mat_traits<M>::read_element_idx(i,j,a);
0052                     ++j2;
0053                     }
0054                 }
0055             T d=determinant_impl_(m,det_size<N-1>());
0056             if( j1&1 )
0057                 d=-d;
0058             det += mat_traits<M>::read_element_idx(0,j1,a) * d;
0059             }
0060         return det;
0061         }
0062
0063     template <class M>
0064     BOOST_QVM_CONSTEXPR BOOST_QVM_INLINE_TRIVIAL
0065     typename mat_traits<M>::scalar_type
0066     determinant_impl( M const & a )
0067         {
0068         BOOST_QVM_STATIC_ASSERT(mat_traits<M>::rows==mat_traits<M>::cols);
0069         return determinant_impl_(a,det_size<mat_traits<M>::rows>());
0070         }
0071     }
0072
0073 } }
0074
0075 #endif