linear_algebra.hpp

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

        >::template apply<
            typename BaseUnitDimensions::item
        > result;
        typedef typename calculate_base_dimension_coefficients_impl<N-1>::template apply<
            typename BaseUnitDimensions::next,
            Dim,
            list<typename result::type, T>
        > next_;
        typedef typename next_::type type;
        typedef list<typename result::next, typename next_::next> next;
    };
};

template<>
struct calculate_base_dimension_coefficients_impl<0> {
    template<class Begin, class BaseUnitDimensions, class T>
    struct apply {
        typedef T type;
        typedef dimensionless_type next;
    };
};

// add_zeroes pushs N zeroes onto the
// front of a list.
//
// list<rational> add_zeroes(list<rational> l, int N) {
//     if(N == 0) {
//         return(l);
//     } else {
//         return(push_front(add_zeroes(l, N-1), 0));
//     }
// }

template<int N>
struct add_zeroes_impl {
    // If you get an error here and your base units are
    // in fact linearly independent, please report it.
    BOOST_MPL_ASSERT_MSG((N > 0), base_units_are_probably_not_linearly_independent, (void));
    template<class T>
    struct apply {
        typedef list<
            static_rational<0>,
            typename add_zeroes_impl<N-1>::template apply<T>::type
        > type;
    };
};

template<>
struct add_zeroes_impl<0> {
    template<class T>
    struct apply {
        typedef T type;
    };
};

// expand_dimensions finds the exponents of
// a set of dimensions in a dimension_list.
// the second parameter is assumed to be
// a superset of the base_dimensions of
// the first parameter.
//
// list<rational> expand_dimensions(dimension_list, list<base_dimension>);

template<int N>
struct expand_dimensions {
    template<class Begin, class DimensionIterator>
    struct apply {
        typedef typename calculate_base_dimension_coefficients_func<
            begins_with_dimension<DimensionIterator>::template apply<typename Begin::item>::value
        >::template apply<DimensionIterator> result;
        typedef list<
            typename result::type,
            typename expand_dimensions<N-1>::template apply<typename Begin::next, typename result::next>::type
        > type;
    };
};

template<>
struct expand_dimensions<0> {
    template<class Begin, class DimensionIterator>
    struct apply {
        typedef dimensionless_type type;
    };
};

template<int N>
struct create_unit_matrix {
    template<class Begin, class Dimensions>
    struct apply {
        typedef typename create_unit_matrix<N - 1>::template apply<typename Begin::next, Dimensions>::type next;
        typedef list<typename expand_dimensions<Dimensions::size::value>::template apply<Dimensions, typename Begin::item::dimension_type>::type, next> type;
    };
};

template<>
struct create_unit_matrix<0> {
    template<class Begin, class Dimensions>
    struct apply {
        typedef dimensionless_type type;
    };
};

template<class T>
struct normalize_units {
    typedef typename find_base_dimensions<T>::type dimensions;
    typedef typename create_unit_matrix<(T::size::value)>::template apply<
        T,
        dimensions
    >::type matrix;
    typedef typename make_square_and_invert<matrix>::type type;
    static const long extra = (type::size::value) - (T::size::value);
};

// multiply_add_units computes M x V
// where M is a matrix and V is a horizontal
// vector
//
// list<rational> multiply_add_units(list<list<rational> >, list<rational>);

template<int N>
struct multiply_add_units_impl {
    template<class Begin1, class Begin2 ,class X>
    struct apply {
        typedef list<
            typename mpl::plus<
                typename mpl::times<
                    typename Begin2::item,
                    X
                >::type,
                typename Begin1::item
            >::type,
            typename multiply_add_units_impl<N-1>::template apply<
                typename Begin1::next,
                typename Begin2::next,
                X
            >::type
        > type;
    };
};

template<>
struct multiply_add_units_impl<0> {
    template<class Begin1, class Begin2 ,class X>
    struct apply {
        typedef dimensionless_type type;
    };
};

template<int N>
struct multiply_add_units {
    template<class Begin1, class Begin2>
    struct apply {
        typedef typename multiply_add_units_impl<
            (Begin2::item::size::value)
        >::template apply<
            typename multiply_add_units<N-1>::template apply<
                typename Begin1::next,
                typename Begin2::next
            >::type,
            typename Begin2::item,
            typename Begin1::item
        >::type type;
    };
};

template<>
struct multiply_add_units<1> {
    template<class Begin1, class Begin2>
    struct apply {
        typedef typename add_zeroes_impl<
            (Begin2::item::size::value)
        >::template apply<dimensionless_type>::type type1;
        typedef typename multiply_add_units_impl<
            (Begin2::item::size::value)
        >::template apply<
            type1,
            typename Begin2::item,
            typename Begin1::item
        >::type type;
    };
};


// strip_zeroes erases the first N elements of a list if
// they are all zero, otherwise returns inconsistent
//
// list strip_zeroes(list l, int N) {
//     if(N == 0) {
//         return(l);
//     } else if(l.front == 0) {
//         return(strip_zeroes(pop_front(l), N-1));
//     } else {
//         return(inconsistent);
//     }
// }

template<int N>
struct strip_zeroes_impl;

template<class T>
struct strip_zeroes_func {
    template<class L, int N>
    struct apply {
        typedef inconsistent type;
    };
};

template<>
struct strip_zeroes_func<static_rational<0> > {
    template<class L, int N>
    struct apply {
        typedef typename strip_zeroes_impl<N-1>::template apply<typename L::next>::type type;
    };
};

template<int N>
struct strip_zeroes_impl {
    template<class T>
    struct apply {
        typedef typename strip_zeroes_func<typename T::item>::template apply<T, N>::type type;
    };
};

template<>
struct strip_zeroes_impl<0> {
    template<class T>
    struct apply {
        typedef T type;
    };
};

// Given a list of base_units, computes the
// exponents of each base unit for a given
// dimension.
//
// list<rational> calculate_base_unit_exponents(list<base_unit> units, dimension_list dimensions);
//
// What is the purpose of all this magic with
// base_dimensions? Can't we just solve the
// equations for the dimension directly?  Yes,
// we can, but remember that solving a
// system of linear equations is O(N^3).
// By normalizing the system we incur a
// high one time cost O(N^4), but for all
// solutions after the first it is O(N^2)
// In addition, the constant factor is
// good because everything is already set up.
// Since we expect a few systems to be
// used many times, the cost of creating
// a system is probably not significant.

template<class T>
struct is_base_dimension_unit {
    typedef mpl::false_ type;
    typedef void base_dimension_type;
};
template<class T>
struct is_base_dimension_unit<list<dim<T, static_rational<1> >, dimensionless_type> > {
    typedef mpl::true_ type;
    typedef T base_dimension_type;
};

template<int N>
struct is_simple_system_impl {
    template<class Begin, class Prev>
    struct apply {
        typedef is_base_dimension_unit<typename Begin::item::dimension_type> test;
        typedef mpl::and_<
            typename test::type,
            mpl::less<Prev, typename test::base_dimension_type>,
            typename is_simple_system_impl<N-1>::template apply<
                typename Begin::next,
                typename test::base_dimension_type
            >
        > type;
        static const bool value = (type::value);
    };
};

template<>
struct is_simple_system_impl<0> {
    template<class Begin, class Prev>
    struct apply : mpl::true_ {
    };
};

template<class T>
struct is_simple_system {
    typedef T Begin;
    typedef is_base_dimension_unit<typename Begin::item::dimension_type> test;
    typedef typename mpl::and_<
        typename test::type,
        typename is_simple_system_impl<
            T::size::value - 1
        >::template apply<
            typename Begin::next::type,
            typename test::base_dimension_type
        >
    >::type type;
    static const bool value = type::value;
};

template<bool>
struct calculate_base_unit_exponents_impl;

template<>
struct calculate_base_unit_exponents_impl<true> {
    template<class T, class Dimensions>
    struct apply {
        typedef typename expand_dimensions<(T::size::value)>::template apply<
            typename find_base_dimensions<T>::type,
            Dimensions
        >::type type;
    };
};

template<>
struct calculate_base_unit_exponents_impl<false> {
    template<class T, class Dimensions>
    struct apply {
        // find the units that correspond to each base dimension
        typedef normalize_units<T> base_solutions;
        // pad the dimension with zeroes so it can just be a
        // list of numbers, making the multiplication easy
        // e.g. if the arguments are list<pound, foot> and
        // list<mass,time^-2> then this step will
        // yield list<0,1,-2>
        typedef typename expand_dimensions<(base_solutions::dimensions::size::value)>::template apply<
            typename base_solutions::dimensions,
            Dimensions
        >::type dimensions;
        // take the unit corresponding to each base unit
        // multiply each of its exponents by the exponent
        // of the base_dimension in the result and sum.
        typedef typename multiply_add_units<dimensions::size::value>::template apply<
            dimensions,
            typename base_solutions::type
        >::type units;
        // Now, verify that the dummy units really
        // cancel out and remove them.
        typedef typename strip_zeroes_impl<base_solutions::extra>::template apply<units>::type type;
    };
};

template<class T, class Dimensions>
struct calculate_base_unit_exponents {
    typedef typename calculate_base_unit_exponents_impl<is_simple_system<T>::value>::template apply<T, Dimensions>::type type;
};

} // namespace detail

} // namespace units

} // namespace boost

#endif