concepts.hpp

矩阵运算源码最新版本

// Copyright 2006. Peter Gottschling, Matthias Troyer, Rolf Bonderer
// Software License for MTL
// 
// Copyright (c) 2007 The Trustees of Indiana University. All rights reserved.
// Authors: Peter Gottschling and Andrew Lumsdaine
// 
// This file is part of the Matrix Template Library
// 
// See also license.mtl.txt in the distribution.

#ifndef LA_CONCEPTS_INCLUDE
#define LA_CONCEPTS_INCLUDE

#include <boost/config/concept_macros.hpp>

#ifdef __GXX_CONCEPTS__
#  include <concepts>
#else
#  ifdef LA_SHOW_WARNINGS
#    warning "Concepts are not used"
#  endif
#endif


#include <boost/numeric/linear_algebra/identity.hpp>
#include <boost/numeric/linear_algebra/is_invertible.hpp>
#include <boost/numeric/linear_algebra/inverse.hpp>
#include <boost/numeric/linear_algebra/operators.hpp>
#include <boost/numeric/linear_algebra/algebraic_concepts.hpp>
#include <complex>

// If desired one can disable the default concept maps with LA_NO_CONCEPT_MAPS

// We consider to change the namespace from math to numeric
// More precisely, the concepts may be moved into namespace numeric and the standard functions stay in math

/// Namespace for mathematical concepts
/** In contrast to the ones in algebra the concepts can require basic implementation concepts like std::Assignable */  
namespace math {

#ifdef __GXX_CONCEPTS__

// ==================================
// Classification of Arithmetic Types
// ==================================

// In addtion to std::Integral
concept Float<typename T> 
  : std::DefaultConstructible<T>, std::CopyConstructible<T>,
    std::LessThanComparable<T>, std::EqualityComparable<T>
{
  T operator+(T);
  T operator+(T, T);
  T& operator+=(T&, T);
  T operator-(T, T);
  T operator-(T);
  T& operator-=(T&, T);
  T operator*(T, T);
  T& operator*=(T&, T);
  T operator/(T, T);
  T& operator/=(T&, T);

  // TBD: Some day, these will come from LessThanComparable,
  // EqualityComparable, etc.
  bool operator>(T, T);
  bool operator<=(T, T);
  bool operator>=(T, T);
  bool operator!=(T, T);

  requires std::Assignable<T>, std::SameType<std::Assignable<T>::result_type, T&>;
}

concept_map Float<float> {}
concept_map Float<double> {}
concept_map Float<long double> {}

// The difference to Float is the lack of LessThanComparable
concept Complex<typename T> 
  : std::DefaultConstructible<T>, std::CopyConstructible<T>,
    std::EqualityComparable<T>
{
  T operator+(T);
  T operator+(T, T);
  T& operator+=(T&, T);
  T operator-(T, T);
  T operator-(T);
  T& operator-=(T&, T);
  T operator*(T, T);
  T& operator*=(T&, T);
  T operator/(T, T);
  T& operator/=(T&, T);

  // TBD: Some day, these will come from EqualityComparable
  bool operator!=(T, T);

  requires std::Assignable<T>, std::SameType<std::Assignable<T>::result_type, T&>;
}

template <typename T>
  requires Float<T>
concept_map Complex<std::complex<T> > {}


// TBD: Concept Arithmetic is useless like this, it should have operations and then be the base for
// Integral, Float and Complex
concept Arithmetic<typename T> {}

template <typename T>
  requires std::Integral<T>
concept_map Arithmetic<T> {}

template <typename T>
  requires Float<T>
concept_map Arithmetic<T> {}

template <typename T>
  requires Arithmetic<T>
concept_map Arithmetic< std::complex<T> > {}




// ================
// Utility Concepts
// ================

// Concepts for functions mapping to same type or convertible
auto concept UnaryIsoFunction<typename Operation, typename Element>
{
    requires std::Callable1<Operation, Element>;
    requires std::Convertible<std::Callable1<Operation, Element>::result_type, Element>;

    typename result_type = std::Callable1<Operation, Element>::result_type;
};


auto concept BinaryIsoFunction<typename Operation, typename Element>
{
    requires std::Callable2<Operation, Element, Element>;
    requires std::Convertible<std::Callable2<Operation, Element, Element>::result_type, Element>;

    typename result_type = std::Callable2<Operation, Element, Element>::result_type;
};

#if 0
auto concept CompatibleBinaryFunction<typename A1, typename A2, typename Result>
{
    typename result_type;
    result_type F(A1, A2);
    requires std::Convertible<result_type, Result>;
}
#endif

// ==================
// Algebraic Concepts
// ==================


auto concept Magma<typename Operation, typename Element>
    : BinaryIsoFunction<Operation, Element>
{
    requires std::Assignable<Element>;
    requires std::Assignable<Element, BinaryIsoFunction<Operation, Element>::result_type>;
};


// For algebraic structures that are commutative but not associative
// As an example floating point numbers are commutative but not associative
//   w.r.t. addition and multiplication
auto concept CommutativeMagma<typename Operation, typename Element>
  : Magma<Operation, Element>, 
    algebra::Commutative<Operation, Element>
{};


// SemiGroup is a refinement which must be nominal
auto concept SemiGroup<typename Operation, typename Element>
  : Magma<Operation, Element>, 
    algebra::SemiGroup<Operation, Element>
{};


auto concept CommutativeSemiGroup<typename Operation, typename Element>
  : SemiGroup<Operation, Element>,
    CommutativeMagma<Operation, Element>
{};


// Adding identity
// auto 
concept Monoid<typename Operation, typename Element>
  : SemiGroup<Operation, Element>, 
    algebra::Monoid<Operation, Element> 
{
    requires std::Convertible<identity_result_type, Element>;
};


auto concept CommutativeMonoid<typename Operation, typename Element>
  : CommutativeSemiGroup<Operation, Element>, 
    Monoid<Operation, Element>
{};


concept PartiallyInvertibleMonoid<typename Operation, typename Element>
  : Monoid<Operation, Element>, 
    algebra::Inversion<Operation, Element> 
{
    typename is_invertible_result_type;
    is_invertible_result_type is_invertible(Operation, Element);
    requires std::Convertible<is_invertible_result_type, bool>;

    requires std::Convertible<inverse_result_type, Element>;

    // Does it overwrites the axiom from algebra::Inversion
    axiom Inversion(Operation op, Element x)
    {
	// Only for invertible elements:
	if (is_invertible(op, x))
	    op( x, inverse(op, x) ) == identity(op, x); 
	if ( is_invertible(op, x) )
	    op( inverse(op, x), x ) == identity(op, x); 
    }
};


auto concept PartiallyInvertibleCommutativeMonoid<typename Operation, typename Element>
  : PartiallyInvertibleMonoid<Operation, Element>, 
    CommutativeMonoid<Operation, Element>   
{};


concept Group<typename Operation, typename Element>
  : PartiallyInvertibleMonoid<Operation, Element>,
    algebra::Group<Operation, Element>
{
    axiom AlwaysInvertible(Operation op, Element x)
    {
	is_invertible(op, x);
    }

    axiom Inversion(Operation op, Element x)
    {
	// In fact this is implied by AlwaysInvertible and inherited Inversion axiom
	// However, we don't rely on the compiler to deduce this
	op( x, inverse(op, x) ) == identity(op, x);
	op( inverse(op, x), x ) == identity(op, x);
    }
};


auto concept AbelianGroup<typename Operation, typename Element>
  : Group<Operation, Element>, 
    PartiallyInvertibleCommutativeMonoid<Operation, Element>,
    algebra::AbelianGroup<Operation, Element>
{};


// ========================
// Additive scalar concepts
// ========================


concept AdditiveMagma<typename Element>
  : Magma< math::add<Element>, Element >
{
    typename plus_assign_result_type;  
    plus_assign_result_type operator+=(Element& x, Element y);
    // requires std::Convertible<plus_assign_result_type, Element>;

    // Operator + is by default defined with +=
    typename addition_result_type;  
    addition_result_type operator+(Element x, Element y);
#if 0
    {
	Element tmp(x);
	return tmp += y;                      defaults NYS
    }
#endif 
    requires std::Convertible<addition_result_type, Element>;

    // Type consistency with Magma
    requires std::Convertible< addition_result_type,   
                            Magma< math::add<Element>, Element >::result_type>;

    // SameType requires more rigorous specializations on pure algebraic functors
    // requires std::SameType< addition_result_type, 
    // 	                    Magma< math::add<Element>, Element >::result_type>;

    axiom Consistency(math::add<Element> op, Element x, Element y)
    {
	op(x, y) == x + y;     
               
	// Consistency definition between + and += might change later
        x + y == x += y;
	// Element tmp = x; tmp+= y; tmp == x + y; not proposal-compliant
    }   
}


auto concept AdditiveCommutativeMagma<typename Element>
  : AdditiveMagma<Element>,
    CommutativeMagma< math::add<Element>, Element >
{};


auto concept AdditiveSemiGroup<typename Element>
  : AdditiveMagma<Element>, 
    SemiGroup< math::add<Element>, Element >
{};


// We really need only one of the additive concepts for the requirements, 
// the requirements of the other would be implied.
// Vice versa, to derive concept maps of nested concepts from
// concept maps of refined concepts, they are needed all.
auto concept AdditiveCommutativeSemiGroup<typename Element>
  : AdditiveSemiGroup<Element>,
    AdditiveCommutativeMagma<Element>,
    CommutativeSemiGroup< math::add<Element>, Element >
{};


concept AdditiveMonoid<typename Element>
  : AdditiveSemiGroup<Element>,
    Monoid< math::add<Element>, Element >
{
    Element zero(Element v);

    axiom Consistency (math::add<Element> op, Element x)
    {
	zero(x) == identity(op, x);
    }
};


// We really need only one of the additive concepts for the requirements, 
// the requirements of the other would be implied.
// Vice versa, to derive concept maps of nested concepts from
// concept maps of refined concepts, they are needed all.
auto concept AdditiveCommutativeMonoid<typename Element>
  : AdditiveMonoid<Element>,
    AdditiveCommutativeSemiGroup<Element>,
    CommutativeMonoid< math::add<Element>, Element >
{};


concept AdditivePartiallyInvertibleMonoid<typename Element>
  : AdditiveMonoid<Element>,
    PartiallyInvertibleMonoid< math::add<Element>, Element >
{
    typename minus_assign_result_type;  
    minus_assign_result_type operator-=(Element& x, Element y);
    // requires std::Convertible<minus_assign_result_type, Element>;
     
    // Operator - by default defined with -=
    typename subtraction_result_type;  
    subtraction_result_type operator-(Element& x, Element y);
#if 0
    {
	Element tmp(x);
	return tmp -= y;                      defaults NYS
    }
#endif 
    requires std::Convertible<subtraction_result_type, Element>;


    typename unary_result_type;  
    unary_result_type operator-(Element x);
#if 0
    {
	return zero(x) - x;      defaults NYS
    }
#endif 
    requires std::Convertible<unary_result_type, Element>;
    
    axiom Consistency(math::add<Element> op, Element x, Element y)
    {
	// consistency between additive and pure algebraic concept
	if ( is_invertible(op, y) )
	    op(x, inverse(op, y)) == x - y;            
	if ( is_invertible(op, y) )
	    inverse(op, y) == -y;                      

	// consistency between unary and binary -
	if ( is_invertible(op, x) )
	    identity(op, x) - x == -x;                 

	// Might change later
	if ( is_invertible(op, y) )
	    x - y == x -= y;                                                       
	// Element tmp = x; tmp-= y; tmp == x - y; not proposal-compliant
    }  

};


auto concept AdditivePartiallyInvertibleCommutativeMonoid<typename Element>
  : AdditivePartiallyInvertibleMonoid<Element>,
    AdditiveCommutativeMonoid<Element>, 
    PartiallyInvertibleCommutativeMonoid< math::add<Element>, Element >
{};



auto concept AdditiveGroup<typename Element>
  : AdditivePartiallyInvertibleMonoid<Element>,
    Group< math::add<Element>, Element >
{};


auto concept AdditiveAbelianGroup<typename Element>
  : AdditiveGroup<Element>,
    AdditiveCommutativeMonoid<Element>,
    AbelianGroup< math::add<Element>, Element >
{};


// ============================
// Multiplitive scalar concepts
// ============================


concept MultiplicativeMagma<typename Element>
  : Magma< math::mult<Element>, Element >
{
    typename mult_assign_result_type;  
    mult_assign_result_type operator*=(Element& x, Element y);
    // requires std::Convertible<mult_assign_result_type, Element>;

    // Operator * is by default defined with *=
    typename mult_result_type;  
    mult_result_type operator*(Element x, Element y);
#if 0
    {
	Element tmp(x);
	return tmp *= y;                      defaults NYS
    }
#endif 
    requires std::Convertible<mult_result_type, Element>;
    
    // Type consistency with Magma
    requires std::Convertible< mult_result_type,   
                            Magma< math::mult<Element>, Element >::result_type>;

    // SameType requires more rigorous specializations on pure algebraic functors
    // requires std::SameType< mult_result_type, 
    // 	                    Magma< math::mult<Element>, Element >::result_type>;


    axiom Consistency(math::mult<Element> op, Element x, Element y)
    {
	op(x, y) == x * y;                 
   
	// Consistency definition between * and *= might change later
        x * y == x *= y;
	// Element tmp = x; tmp*= y; tmp == x * y; not proposal-compliant
    }  

}


auto concept MultiplicativeSemiGroup<typename Element>
  : MultiplicativeMagma<Element>,
    SemiGroup< math::mult<Element>, Element >
{};


auto concept MultiplicativeCommutativeSemiGroup<typename Element>
  : MultiplicativeSemiGroup<Element>,
    CommutativeSemiGroup< math::mult<Element>, Element >
{};


concept MultiplicativeMonoid<typename Element>
  : MultiplicativeSemiGroup<Element>,
    Monoid< math::mult<Element>, Element >
{
    Element one(Element v);

    axiom Consistency (math::mult<Element> op, Element x)
    {
	one(x) == identity(op, x);
    }
};


auto concept MultiplicativeCommutativeMonoid<typename Element>
  : MultiplicativeMonoid<Element>,
    MultiplicativeCommutativeSemiGroup<Element>,
    CommutativeMonoid< math::mult<Element>, Element >
{};


concept MultiplicativePartiallyInvertibleMonoid<typename Element>
  : MultiplicativeMonoid<Element>,
    PartiallyInvertibleMonoid< math::mult<Element>, Element >
{
    typename divide_assign_result_type;  
    divide_assign_result_type operator/=(Element& x, Element y);
    // requires std::Convertible<divide_assign_result_type, Element>;
     
    // Operator / by default defined with /=
    typename division_result_type = Element;  
    division_result_type operator/(Element x, Element y);
#if 0
    {
	Element tmp(x);
	return tmp /= y;                      defaults NYS
    }
#endif 
    requires std::Convertible<division_result_type, Element>;
    
    axiom Consistency(math::mult<Element> op, Element x, Element y)
    {
	// consistency between multiplicative and pure algebraic concept
	if ( is_invertible(op, y) )
	    op(x, inverse(op, y)) == x / y;            

	// Consistency between / and /=, might change later
	if ( is_invertible(op, y) )
	    x / y == x /= y;              
	// Element tmp = x; tmp/= y; tmp == x / y; not proposal-compliant 
    }  
};
 

auto concept MultiplicativePartiallyInvertibleCommutativeMonoid<typename Element>
  : MultiplicativePartiallyInvertibleMonoid<Element>,
    MultiplicativeCommutativeMonoid<Element>,
    PartiallyInvertibleCommutativeMonoid< math::mult<Element>, Element >
{};
 

auto concept MultiplicativeGroup<typename Element>
  : MultiplicativeMonoid<Element>,
    Group< math::mult<Element>, Element >
{};


auto concept MultiplicativeAbelianGroup<typename Element>
  : MultiplicativeGroup<Element>,
    MultiplicativeCommutativeMonoid<Element>,
    AbelianGroup< math::mult<Element>, Element >
{};