linear_algebra.hpp

来自「Boost provides free peer-reviewed portab」· HPP 代码 · 共 1,075 行 · 第 1/3 页

    struct apply {
        typedef typename RowsBegin::item current_row;
        typedef typename current_row::item current_value;
        typedef typename divide_equation<(current_row::size::value - 1)>::template apply<typename current_row::next, current_value>::type new_equation;
        typedef typename divide_equation<(IdentityBegin::item::size::value)>::template apply<typename IdentityBegin::item, current_value>::type transformed_identity_equation;
        typedef typename invert_handle_after_pivot_row<(RowsBegin::size::value - 1)>::template apply<
            typename RowsBegin::next,
            typename IdentityBegin::next,
            new_equation,
            transformed_identity_equation
        > next;

        // results
        // Note that we don't add the pivot row to the
        // results here, because it needs to propagated up
        // to the diagonal.
        typedef typename next::new_matrix new_matrix;
        typedef typename next::identity_result identity_result;
        typedef new_equation pivot_row;
        typedef transformed_identity_equation identity_pivot_row;
    };
};

// The one and only non-zero element--at the end
template<>
struct invert_strip_leading_zeroes<false, true> {
    template<class RowsBegin, class IdentityBegin>
    struct apply {
        typedef typename RowsBegin::item current_row;
        typedef typename current_row::item current_value;
        typedef typename divide_equation<(current_row::size::value - 1)>::template apply<typename current_row::next, current_value>::type new_equation;
        typedef typename divide_equation<(IdentityBegin::item::size::value)>::template apply<typename IdentityBegin::item, current_value>::type transformed_identity_equation;

        // results
        // Note that we don't add the pivot row to the
        // results here, because it needs to propagated up
        // to the diagonal.
        typedef dimensionless_type identity_result;
        typedef dimensionless_type new_matrix;
        typedef new_equation pivot_row;
        typedef transformed_identity_equation identity_pivot_row;
    };
};

// One of the initial zeroes
template<>
struct invert_strip_leading_zeroes<true, false> {
    template<class RowsBegin, class IdentityBegin>
    struct apply {
        typedef typename RowsBegin::item current_row;
        typedef typename RowsBegin::next::item next_row;
        typedef typename invert_strip_leading_zeroes<
            next_row::item::Numerator == 0,
            RowsBegin::size::value == 2
        >::template apply<
            typename RowsBegin::next,
            typename IdentityBegin::next
        > next;
        typedef typename IdentityBegin::item current_identity_row;
        // these are propagated up.
        typedef typename next::pivot_row pivot_row;
        typedef typename next::identity_pivot_row identity_pivot_row;
        typedef list<
            typename eliminate_from_pair_of_equations_impl<(current_row::size::value - 1)>::template apply<
                typename current_row::next,
                pivot_row,
                typename current_row::item,
                static_rational<1>
            >::type,
            typename next::new_matrix
        > new_matrix;
        typedef list<
            typename eliminate_from_pair_of_equations_impl<(current_identity_row::size::value)>::template apply<
                current_identity_row,
                identity_pivot_row,
                typename current_row::item,
                static_rational<1>
            >::type,
            typename next::identity_result
        > identity_result;
    };
};

// the last element, and is zero.
// Should never happen.
template<>
struct invert_strip_leading_zeroes<true, true> {
};

template<int N>
struct invert_handle_after_pivot_row {
    template<class RowsBegin, class IdentityBegin, class MatrixPivot, class IdentityPivot>
    struct apply {
        typedef typename invert_handle_after_pivot_row<N - 1>::template apply<
            typename RowsBegin::next,
            typename IdentityBegin::next,
            MatrixPivot,
            IdentityPivot
        > next;
        typedef typename RowsBegin::item current_row;
        typedef typename IdentityBegin::item current_identity_row;
        typedef MatrixPivot pivot_row;
        typedef IdentityPivot identity_pivot_row;

        // results
        typedef list<
            typename eliminate_from_pair_of_equations_impl<(current_row::size::value - 1)>::template apply<
                typename current_row::next,
                pivot_row,
                typename current_row::item,
                static_rational<1>
            >::type,
            typename next::new_matrix
        > new_matrix;
        typedef list<
            typename eliminate_from_pair_of_equations_impl<(current_identity_row::size::value)>::template apply<
                current_identity_row,
                identity_pivot_row,
                typename current_row::item,
                static_rational<1>
            >::type,
            typename next::identity_result
        > identity_result;
    };
};

template<>
struct invert_handle_after_pivot_row<0> {
    template<class RowsBegin, class IdentityBegin, class MatrixPivot, class IdentityPivot>
    struct apply {
        typedef dimensionless_type new_matrix;
        typedef dimensionless_type identity_result;
    };
};

template<int N>
struct invert_impl {
    template<class RowsBegin, class IdentityBegin>
    struct apply {
        typedef typename invert_handle_inital_rows<RowsBegin::size::value - N>::template apply<RowsBegin, IdentityBegin> process_column;
        typedef typename invert_impl<N - 1>::template apply<
            typename process_column::new_matrix,
            typename process_column::identity_result
        >::type type;
    };
};

template<>
struct invert_impl<0> {
    template<class RowsBegin, class IdentityBegin>
    struct apply {
        typedef IdentityBegin type;
    };
};

template<int N>
struct make_identity {
    template<int Size>
    struct apply {
        typedef list<typename create_row_of_identity<Size - N, Size>::type, typename make_identity<N - 1>::template apply<Size>::type> type;
    };
};

template<>
struct make_identity<0> {
    template<int Size>
    struct apply {
        typedef dimensionless_type type;
    };
};

template<class Matrix>
struct make_square_and_invert {
    typedef typename Matrix::item top_row;
    typedef typename determine_extra_equations<(top_row::size::value), false>::template apply<
        Matrix,                 // RowsBegin
        top_row::size::value,   // TotalColumns
        Matrix                  // Result
    >::type invertible;
    typedef typename invert_impl<invertible::size::value>::template apply<
        invertible,
        typename make_identity<invertible::size::value>::template apply<invertible::size::value>::type
    >::type type;
};


// find_base_dimensions takes a list of
// base_units and returns a sorted list
// of all the base_dimensions they use.
//
// list<base_dimension> find_base_dimensions(list<base_unit> l) {
//     set<base_dimension> dimensions;
//     for_each(base_unit unit : l) {
//         for_each(dim d : unit.dimension_type) {
//             dimensions = insert(dimensions, d.tag_type);
//         }
//     }
//     return(sort(dimensions, _1 > _2, front_inserter(list<base_dimension>())));
// }

typedef char set_no;
struct set_yes { set_no dummy[2]; };

template<class T>
struct wrap {};

struct set_end {
    static set_no lookup(...);
    typedef mpl::long_<0> size;
};

template<class T, class Next>
struct set : Next {
    using Next::lookup;
    static set_yes lookup(wrap<T>*);
    typedef T item;
    typedef Next next;
    typedef typename mpl::next<typename Next::size>::type size;
};

template<bool has_key>
struct set_insert;

template<>
struct set_insert<true> {
    template<class Set, class T>
    struct apply {
        typedef Set type;
    };
};

template<>
struct set_insert<false> {
    template<class Set, class T>
    struct apply {
        typedef set<T, Set> type;
    };
};

template<class Set, class T>
struct has_key {
    static const long size = sizeof(Set::lookup((wrap<T>*)0));
    static const bool value = (size == sizeof(set_yes));
};

template<int N>
struct find_base_dimensions_impl_impl {
    template<class Begin, class S>
    struct apply {
        typedef typename find_base_dimensions_impl_impl<N-1>::template apply<
            typename Begin::next,
            S
        >::type next;

        typedef typename set_insert<
            (has_key<next, typename Begin::item::tag_type>::value)
        >::template apply<
            next,
            typename Begin::item::tag_type
        >::type type;
    };
};

template<>
struct find_base_dimensions_impl_impl<0> {
    template<class Begin, class S>
    struct apply {
        typedef S type;
    };
};

template<int N>
struct find_base_dimensions_impl {
    template<class Begin>
    struct apply {
        typedef typename find_base_dimensions_impl_impl<(Begin::item::dimension_type::size::value)>::template apply<
            typename Begin::item::dimension_type,
            typename find_base_dimensions_impl<N-1>::template apply<typename Begin::next>::type
        >::type type;
    };
};

template<>
struct find_base_dimensions_impl<0> {
    template<class Begin>
    struct apply {
        typedef set_end type;
    };
};

template<class T>
struct find_base_dimensions {
    typedef typename insertion_sort<
        typename find_base_dimensions_impl<
            (T::size::value)
        >::template apply<T>::type
    >::type type;
};

// calculate_base_dimension_coefficients finds
// the coefficients corresponding to the first
// base_dimension in each of the dimension_lists.
// It returns two values.  The first result
// is a list of the coefficients.  The second
// is a list with all the incremented iterators.
// When we encounter a base_dimension that is
// missing from a dimension_list, we do not
// increment the iterator and we set the
// coefficient to zero.

template<bool has_dimension>
struct calculate_base_dimension_coefficients_func;

template<>
struct calculate_base_dimension_coefficients_func<true> {
    template<class T>
    struct apply {
        typedef typename T::item::value_type type;
        typedef typename T::next next;
    };
};

template<>
struct calculate_base_dimension_coefficients_func<false> {
    template<class T>
    struct apply {
        typedef static_rational<0> type;
        typedef T next;
    };
};

// begins_with_dimension returns true iff its first
// parameter is a valid iterator which yields its
// second parameter when dereferenced.

template<class Iterator>
struct begins_with_dimension {
    template<class Dim>
    struct apply : 
        boost::is_same<
            Dim,
            typename Iterator::item::tag_type
        > {};
};

template<>
struct begins_with_dimension<dimensionless_type> {
    template<class Dim>
    struct apply : mpl::false_ {};
};

template<int N>
struct calculate_base_dimension_coefficients_impl {
    template<class BaseUnitDimensions,class Dim,class T>
    struct apply {
        typedef typename calculate_base_dimension_coefficients_func<
            begins_with_dimension<typename BaseUnitDimensions::item>::template apply<
                Dim
            >::value