LCPQ/dft_tools
3
0
Fork 0
mirror of https://github.com/triqs/dft_tools synced 2024-12-25 05:43:40 +01:00
Code Issues Releases Wiki Activity

arrays : revised version of map & fold

Browse Source 
- also cleaned a bit reinterpretview example
This commit is contained in:
Olivier Parcollet 2013-09-10 11:18:04 +02:00
parent a0eba01769
commit c06e183ae5
12 changed files with 274 additions and 285 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
test/triqs/arrays/fold_expression.cpp
View File

@ -28,7 +28,7 @@
using namespace triqs::arrays; using namespace triqs::arrays;
auto Abs = map( boost::function< double (double)> ( static_cast< double (*)(double)> (std::abs)) ); auto Abs = map( std::function< double (double)> ( static_cast< double (*)(double)> (std::abs)) );
int main(int argc, char **argv) { int main(int argc, char **argv) {

9
test/triqs/arrays/group_indices.cpp
View File

@ -5,7 +5,6 @@
namespace tqa=triqs::arrays; namespace tqa=triqs::arrays;
namespace tql=triqs::clef; namespace tql=triqs::clef;
namespace mpl=boost::mpl;
using tqa::m_index; using tqa::m_index;
@ -57,7 +56,9 @@ int main () {
test<TRAVERSAL_ORDER_FORTRAN>(); test<TRAVERSAL_ORDER_FORTRAN>();
test<TRAVERSAL_ORDER_C>(); test<TRAVERSAL_ORDER_C>();
// an alternative syntax ? Why to use the npl here ?? // template dans le namespace clef pour ADL et enable sur le concept de A ...
//auto V = tqa::group_indices( A(i_,j_,k_,l_), m_index(i_,j_), m_index(k_,l_)); // an alternative syntax ?
//auto V = group_indices( A(i_,j_,k_,l_), var(i_,j_), var(k_,l_));
// match A(...), deduce the int of hte ph, and make a finder for a ph in it
// finder<Expr>::find<i>::value --> m_index<...>
} }

18
test/triqs/arrays/mapped_functions.cpp
View File

@ -22,13 +22,14 @@
#include "./src/mapped_functions.hpp" #include "./src/mapped_functions.hpp"
#include "./src/matrix.hpp" #include "./src/matrix.hpp"
#include "./src/array.hpp"
#include <iostream> #include <iostream>
namespace tqa = triqs::arrays; using namespace triqs::arrays;
template<typename T> void test( T val=1 ) { template<typename T> void test( T val=1 ) {
tqa::matrix<T> A(3,3, FORTRAN_LAYOUT),B(3,3, FORTRAN_LAYOUT); matrix<T> A(3,3, FORTRAN_LAYOUT),B(3,3, FORTRAN_LAYOUT);
typedef array_view<T,2> V;
for (int i =0; i<3; ++i) for (int i =0; i<3; ++i)
for (int j=0; j<3; ++j) for (int j=0; j<3; ++j)
{ A(i,j) = (i+2*j+1); B(i,j) = (i-j);} { A(i,j) = (i+2*j+1); B(i,j) = (i-j);}
@ -36,8 +37,9 @@ template<typename T> void test( T val=1 ) {
A *=val; B*=val; A *=val; B*=val;
T s = 10; T s = 10;
TEST(A); TEST(A);
TEST(make_matrix(pow(A,2)));
TEST(make_matrix(cosh(A))); TEST(make_matrix(pow(V(A),2)));
TEST(make_matrix(cosh(V(A))));
TEST(B); TEST(B);
TEST(abs(B)); TEST(abs(B));
TEST(make_matrix(abs(B))); TEST(make_matrix(abs(B)));
@ -48,6 +50,12 @@ template<typename T> void test( T val=1 ) {
TEST(make_matrix(real(B))); TEST(make_matrix(real(B)));
TEST(make_matrix(imag(B))); TEST(make_matrix(imag(B)));
auto aa = array<T,2>{ { 1,2}, {3,4}};
TEST(make_matrix(exp(aa)));
// does not compile, since exp is only element wise at the moment
// to do : implement it for matrix...
//TEST(make_matrix(exp( matrix<double>{{ 1,2}, {3,4}} )));
} }
int main(int argc, char **argv) { int main(int argc, char **argv) {

197
test/triqs/arrays/mapped_functions.output
View File

@ -2,93 +2,11 @@
[[1,3,5] [[1,3,5]
[2,4,6] [2,4,6]
[3,5,7]] [3,5,7]]
(make_matrix(pow(A,2))) ---> (make_matrix(pow(V(A),2))) --->
[[1,9,25] [[1,9,25]
[4,16,36] [4,16,36]
[9,25,49]] [9,25,49]]
(make_matrix(cosh(A))) ---> (make_matrix(cosh(V(A)))) --->
[[1,10,74]
[3,27,201]
[10,74,548]]
(B) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(abs(B)) ---> lazy matrix resulting of a mapping
(make_matrix(abs(B))) --->
[[0,1,2]
[1,0,1]
[2,1,0]]
(make_matrix(abs(B+B))) --->
[[0,2,4]
[2,0,2]
[4,2,0]]
(make_matrix(A+ s*B)) --->
[[1,-7,-15]
[12,4,-4]
[23,15,7]]
(make_matrix(abs(A+s*B))) --->
[[1,7,15]
[12,4,4]
[23,15,7]]
(make_matrix(real(B))) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(make_matrix(imag(B))) --->
[[0,0,0]
[0,0,0]
[0,0,0]]
(A) --->
[[1,3,5]
[2,4,6]
[3,5,7]]
(make_matrix(pow(A,2))) --->
[[1,9,25]
[4,16,36]
[9,25,49]]
(make_matrix(cosh(A))) --->
[[1,10,74]
[3,27,201]
[10,74,548]]
(B) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(abs(B)) ---> lazy matrix resulting of a mapping
(make_matrix(abs(B))) --->
[[0,1,2]
[1,0,1]
[2,1,0]]
(make_matrix(abs(B+B))) --->
[[0,2,4]
[2,0,2]
[4,2,0]]
(make_matrix(A+ s*B)) --->
[[1,-7,-15]
[12,4,-4]
[23,15,7]]
(make_matrix(abs(A+s*B))) --->
[[1,7,15]
[12,4,4]
[23,15,7]]
(make_matrix(real(B))) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(make_matrix(imag(B))) --->
[[0,0,0]
[0,0,0]
[0,0,0]]
(A) --->
[[1,3,5]
[2,4,6]
[3,5,7]]
(make_matrix(pow(A,2))) --->
[[1,9,25]
[4,16,36]
[9,25,49]]
(make_matrix(cosh(A))) --->
[[1.54308,10.0677,74.2099] [[1.54308,10.0677,74.2099]
[3.7622,27.3082,201.716] [3.7622,27.3082,201.716]
[10.0677,74.2099,548.317]] [10.0677,74.2099,548.317]]
@ -96,7 +14,7 @@
[[0,-1,-2] [[0,-1,-2]
[1,0,-1] [1,0,-1]
[2,1,0]] [2,1,0]]
(abs(B)) ---> lazy matrix resulting of a mapping (abs(B)) ---> mapping result
(make_matrix(abs(B))) ---> (make_matrix(abs(B))) --->
[[0,1,2] [[0,1,2]
[1,0,1] [1,0,1]
@ -121,15 +39,106 @@
[[0,0,0] [[0,0,0]
[0,0,0] [0,0,0]
[0,0,0]] [0,0,0]]
(make_matrix(exp(aa))) --->
[[2.71828,7.38906]
[20.0855,54.5982]]
(A) --->
[[1,3,5]
[2,4,6]
[3,5,7]]
(make_matrix(pow(V(A),2))) --->
[[1,9,25]
[4,16,36]
[9,25,49]]
(make_matrix(cosh(V(A)))) --->
[[1.54308,10.0677,74.2099]
[3.7622,27.3082,201.716]
[10.0677,74.2099,548.317]]
(B) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(abs(B)) ---> mapping result
(make_matrix(abs(B))) --->
[[0,1,2]
[1,0,1]
[2,1,0]]
(make_matrix(abs(B+B))) --->
[[0,2,4]
[2,0,2]
[4,2,0]]
(make_matrix(A+ s*B)) --->
[[1,-7,-15]
[12,4,-4]
[23,15,7]]
(make_matrix(abs(A+s*B))) --->
[[1,7,15]
[12,4,4]
[23,15,7]]
(make_matrix(real(B))) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(make_matrix(imag(B))) --->
[[0,0,0]
[0,0,0]
[0,0,0]]
(make_matrix(exp(aa))) --->
[[2.71828,7.38906]
[20.0855,54.5982]]
(A) --->
[[1,3,5]
[2,4,6]
[3,5,7]]
(make_matrix(pow(V(A),2))) --->
[[1,9,25]
[4,16,36]
[9,25,49]]
(make_matrix(cosh(V(A)))) --->
[[1.54308,10.0677,74.2099]
[3.7622,27.3082,201.716]
[10.0677,74.2099,548.317]]
(B) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(abs(B)) ---> mapping result
(make_matrix(abs(B))) --->
[[0,1,2]
[1,0,1]
[2,1,0]]
(make_matrix(abs(B+B))) --->
[[0,2,4]
[2,0,2]
[4,2,0]]
(make_matrix(A+ s*B)) --->
[[1,-7,-15]
[12,4,-4]
[23,15,7]]
(make_matrix(abs(A+s*B))) --->
[[1,7,15]
[12,4,4]
[23,15,7]]
(make_matrix(real(B))) --->
[[0,-1,-2]
[1,0,-1]
[2,1,0]]
(make_matrix(imag(B))) --->
[[0,0,0]
[0,0,0]
[0,0,0]]
(make_matrix(exp(aa))) --->
[[2.71828,7.38906]
[20.0855,54.5982]]
(A) ---> (A) --->
[[(1,0),(3,0),(5,0)] [[(1,0),(3,0),(5,0)]
[(2,0),(4,0),(6,0)] [(2,0),(4,0),(6,0)]
[(3,0),(5,0),(7,0)]] [(3,0),(5,0),(7,0)]]
(make_matrix(pow(A,2))) ---> (make_matrix(pow(V(A),2))) --->
[[(1,0),(9,0),(25,0)] [[(1,0),(9,0),(25,0)]
[(4,0),(16,0),(36,0)] [(4,0),(16,0),(36,0)]
[(9,0),(25,0),(49,0)]] [(9,0),(25,0),(49,0)]]
(make_matrix(cosh(A))) ---> (make_matrix(cosh(V(A)))) --->
[[(1.54308,0),(10.0677,0),(74.2099,0)] [[(1.54308,0),(10.0677,0),(74.2099,0)]
[(3.7622,0),(27.3082,0),(201.716,0)] [(3.7622,0),(27.3082,0),(201.716,0)]
[(10.0677,0),(74.2099,0),(548.317,0)]] [(10.0677,0),(74.2099,0),(548.317,0)]]
@ -137,7 +146,7 @@
[[(0,0),(-1,0),(-2,0)] [[(0,0),(-1,0),(-2,0)]
[(1,0),(0,0),(-1,0)] [(1,0),(0,0),(-1,0)]
[(2,0),(1,0),(0,0)]] [(2,0),(1,0),(0,0)]]
(abs(B)) ---> lazy matrix resulting of a mapping (abs(B)) ---> mapping result
(make_matrix(abs(B))) ---> (make_matrix(abs(B))) --->
[[0,1,2] [[0,1,2]
[1,0,1] [1,0,1]
@ -162,15 +171,18 @@
[[0,0,0] [[0,0,0]
[0,0,0] [0,0,0]
[0,0,0]] [0,0,0]]
(make_matrix(exp(aa))) --->
[[(2.71828,0),(7.38906,0)]
[(20.0855,0),(54.5982,0)]]
(A) ---> (A) --->
[[(-2,3),(-6,9),(-10,15)] [[(-2,3),(-6,9),(-10,15)]
[(-4,6),(-8,12),(-12,18)] [(-4,6),(-8,12),(-12,18)]
[(-6,9),(-10,15),(-14,21)]] [(-6,9),(-10,15),(-14,21)]]
(make_matrix(pow(A,2))) ---> (make_matrix(pow(V(A),2))) --->
[[(-5,-12),(-45,-108),(-125,-300)] [[(-5,-12),(-45,-108),(-125,-300)]
[(-20,-48),(-80,-192),(-180,-432)] [(-20,-48),(-80,-192),(-180,-432)]
[(-45,-108),(-125,-300),(-245,-588)]] [(-45,-108),(-125,-300),(-245,-588)]]
(make_matrix(cosh(A))) ---> (make_matrix(cosh(V(A)))) --->
[[(-3.72455,-0.511823),(-183.789,-83.1297),(-8366.62,-7161.77)] [[(-3.72455,-0.511823),(-183.789,-83.1297),(-8366.62,-7161.77)]
[(26.2206,7.62523),(1257.75,799.751),(53734.9,61113.4)] [(26.2206,7.62523),(1257.75,799.751),(53734.9,61113.4)]
[(-183.789,-83.1297),(-8366.62,-7161.77),(-329351,-503083)]] [(-183.789,-83.1297),(-8366.62,-7161.77),(-329351,-503083)]]
@ -178,7 +190,7 @@
[[(0,0),(2,-3),(4,-6)] [[(0,0),(2,-3),(4,-6)]
[(-2,3),(0,0),(2,-3)] [(-2,3),(0,0),(2,-3)]
[(-4,6),(-2,3),(0,0)]] [(-4,6),(-2,3),(0,0)]]
(abs(B)) ---> lazy matrix resulting of a mapping (abs(B)) ---> mapping result
(make_matrix(abs(B))) ---> (make_matrix(abs(B))) --->
[[0,3.60555,7.2111] [[0,3.60555,7.2111]
[3.60555,0,3.60555] [3.60555,0,3.60555]
@ -203,3 +215,6 @@
[[0,-3,-6] [[0,-3,-6]
[3,0,-3] [3,0,-3]
[6,3,0]] [6,3,0]]
(make_matrix(exp(aa))) --->
[[(2.71828,0),(7.38906,0)]
[(20.0855,0),(54.5982,0)]]

40
test/triqs/arrays/reinterpret_view.cpp
View File

@ -2,7 +2,7 @@
* *
* TRIQS: a Toolbox for Research in Interacting Quantum Systems * TRIQS: a Toolbox for Research in Interacting Quantum Systems
* *
* Copyright (C) 2011 by O. Parcollet * Copyright (C) 2013 by O. Parcollet
* *
* TRIQS is free software: you can redistribute it and/or modify it under the * TRIQS is free software: you can redistribute it and/or modify it under the
* terms of the GNU General Public License as published by the Free Software * terms of the GNU General Public License as published by the Free Software
@ -24,18 +24,46 @@
namespace triqs { namespace arrays { namespace triqs { namespace arrays {
template <typename V, int R, ull_t Opt, ull_t To, typename ... I> /* template <typename V, int R, ull_t Opt, ull_t To, typename ... I>
array_view<V, sizeof...(I), Opt> reinterpret_array_view (array<V,R,Opt,To> const & a, I ... index) { array_view<V, sizeof...(I), Opt> reinterpret_array_view (array<V,R,Opt,To> const & a, I ... index) {
static int constexpr rank = sizeof...(I); //static int constexpr rank = sizeof...(I);
typedef array_view<V, rank, Opt, indexmaps::mem_layout::c_order(rank)> return_t; typedef array_view<V, sizeof...(I), Opt, indexmaps::mem_layout::c_order(rank)> return_t;
return return_t (typename return_t::indexmap_type (mini_vector<size_t,rank>(index...)) , a.storage()); //typedef array_view<V, rank, Opt, indexmaps::mem_layout::c_order(rank)> return_t;
return return_t { make_shape(index...), a.storage() };
//return return_t (typename return_t::indexmap_type (mini_vector<size_t,rank>(index...)) , a.storage());
} }
*/
template <typename V, int R, ull_t Opt, ull_t To, typename ... I>
array_view<V, sizeof...(I), Opt> reinterpret (array<V,R,Opt,To> const & a, I ... index) {
return { {make_shape(index...)}, a.storage() };
}
// wrong for views
template <typename V, int R, ull_t Opt, ull_t To, bool B, typename ... I>
array_view<V, sizeof...(I), Opt,indexmaps::mem_layout::c_order(sizeof...(I))> reinterpret_array_view (array_view<V,R,Opt,To,B> const & a, I ... index) {
if (!has_contiguous_data(a)) TRIQS_RUNTIME_ERROR << "reinterpretation failure : data of the view are not contiguous";
return { {make_shape(index...)}, a.storage() };
}
/*
template <typename V, int R, ull_t Opt, ull_t To, typename ... I>
array_view<V, sizeof...(I), Opt,indexmaps::mem_layout::c_order(sizeof...(I))> reinterpret_array_view_c (matrix<V,Opt,To> const & a, I ... index) {
return { {make_shape(index...)}, a.storage() };
}
template <typename V, int R, ull_t Opt, ull_t To, typename ... I>
array_view<V, sizeof...(I), Opt,indexmaps::mem_layout::c_order(sizeof...(I))> reinterpret_array_view_c (vector<V,Opt> const & a, I ... index) {
return { {make_shape(index...)}, a.storage() };
}
*/
}} }}
using namespace triqs::arrays;
int main(int argc, char **argv) { int main(int argc, char **argv) {
triqs::arrays::array<long,1> A = {1,2,3,4,5,6}; triqs::arrays::array<long,1> A = {1,2,3,4,5,6};
std::cout << reinterpret_array_view(A, 2,3); std::cout << reinterpret(A, 2,3)<< std::endl;
//std::cout << reinterpret_array_view(A(range(1,5)), 2,3)<< std::endl;
return 0; return 0;
} }

8
triqs/arrays/algorithms.hpp
View File

@ -39,14 +39,14 @@ namespace triqs { namespace arrays {
template<class A> template<class A>
typename A::value_type max_element(A const &a) { typename A::value_type max_element(A const &a) {
//typedef typename A::value_type const & T; return fold ( static_cast<T(*)(T,T)>(std::max<T>)) ( a, get_first_element(a)); typedef typename A::value_type T;
return fold (std::max<typename A::value_type const & >) (a, get_first_element(a)); return fold ([](T const & a, T const & b) { return std::max(a,b);} ) (a, get_first_element(a));
} }
template<class A> template<class A>
typename A::value_type min_element(A const &a) { typename A::value_type min_element(A const &a) {
// typedef typename A::value_type const & T; return fold ( static_cast<T(*)(T,T)>(std::min<T>)) ( a, get_first_element(a)); typedef typename A::value_type T;
return fold (std::min<typename A::value_type const & >) (a, get_first_element(a)); return fold ([](T const & a, T const & b) { return std::min(a,b);} ) (a, get_first_element(a));
} }
template <class A> template <class A>

2
triqs/arrays/asserts.hpp
View File

@ -30,7 +30,7 @@ namespace triqs { namespace arrays {
template<class ArrayType1, class ArrayType2 > template<class ArrayType1, class ArrayType2 >
void assert_all_close( ArrayType1 const & A, ArrayType2 const & B, double precision, bool relative = false) { void assert_all_close( ArrayType1 const & A, ArrayType2 const & B, double precision, bool relative = false) {
typedef typename ArrayType1::value_type F; typedef typename ArrayType1::value_type F;
auto Abs = map( boost::function<double(F)> (assert_abs<F>) ); auto Abs = map( std::function<double(F)> (assert_abs<F>) );
auto r = max_element (Abs(A-B)); auto r = max_element (Abs(A-B));
auto r2 = max_element (Abs(A) + Abs(B)); auto r2 = max_element (Abs(A) + Abs(B));
if ( r > (relative ? precision * r2 : precision) ) if ( r > (relative ? precision * r2 : precision) )

4
triqs/arrays/blas_lapack/getri.hpp
View File

@ -42,8 +42,8 @@ namespace triqs { namespace arrays { namespace lapack {
} }
} }
inline size_t r_round(double x) { return round(x);} inline size_t r_round(double x) { return std::round(x);}
inline size_t r_round(std::complex<double> x) { return round(real(x));} inline size_t r_round(std::complex<double> x) { return std::round(std::real(x));}
/** /**
* Calls getri on a matrix or view * Calls getri on a matrix or view
* Takes care of making temporary copies if necessary * Takes care of making temporary copies if necessary

35
triqs/arrays/functional/fold.hpp
View File

@ -28,37 +28,24 @@
namespace triqs { namespace arrays { namespace triqs { namespace arrays {
template<class F> template<class F>
class fold_worker { struct fold_worker {
F f; F f;
typedef typename boost::remove_const<typename boost::remove_reference<typename F::result_type >::type>::type result_type;
template<class A> template<class A, class R> struct fold_func_adaptor {
struct fold_func_adaptor { F const & f; A const & a; R & r;
A const & b; template<typename ... Args> void operator()(Args const & ... args) { r = f(r,a(args...)); }
F f; result_type r;
fold_func_adaptor(F f_, A const & a_,typename A::value_type init ):b(a_), f(f_) { r= init;}
template<typename ... Args> void operator()(Args const & ... args) { r = f(r,b(args...)); }
}; };
public: template<class A, class R>
R operator() (A const & a, R init) const {
fold_worker ( F const & f_):f(f_) {} foreach(a, fold_func_adaptor<A,R> {f,a,init});
return init;
template<class A>
result_type operator() (A const & a, typename A::value_type init = typename A::value_type() ) const {
fold_func_adaptor<A> func(f,a,init);
//fold_func_adaptor<typename A::value_type> func(f,init);
foreach(a,std::ref(func));
//foreach(a,boost::ref(func));
return func.r;
} }
template<class A> typename A::value_type operator() (A const & a) const { return (*this)(a, typename A::value_type{});}
}; };
template<class F> fold_worker<F> fold (F const & f) { return fold_worker<F>(f);} template<class F> fold_worker<F> fold (F f) { return {std::move(f)};}
template<typename R, typename A1>
fold_worker< boost::function<R(A1,A1)> >
fold (R (*f)(A1,A1)) { return fold( boost::function<R(A1,A1)>(f)); }
}}//namespace triqs::arrays }}//namespace triqs::arrays

149
triqs/arrays/functional/map.hpp
View File

@ -22,6 +22,8 @@
#define TRIQS_ARRAYS_EXPRESSION_MAP_H #define TRIQS_ARRAYS_EXPRESSION_MAP_H
#include "../impl/common.hpp" #include "../impl/common.hpp"
#include <functional> #include <functional>
//#include "../../utility/function_arg_ret_type.hpp"
namespace triqs { namespace arrays { namespace triqs { namespace arrays {
template<class F, int arity=F::arity> class map_impl; template<class F, int arity=F::arity> class map_impl;
@ -30,105 +32,80 @@ namespace triqs { namespace arrays {
* Given a function f : arg_type -> result_type, map(f) is the function promoted to arrays * Given a function f : arg_type -> result_type, map(f) is the function promoted to arrays
* map(f) : array<arg_type, N, Opt> --> array<result_type, N, Opt> * map(f) : array<arg_type, N, Opt> --> array<result_type, N, Opt>
*/ */
template<class F> map_impl<F,1> map (F const & f) { return map_impl<F,1>(f,true); } //template<class F> map_impl<F,utility::function_arg_ret_type<F>::arity> map (F f) { return {std::move(f),true}; }
template<class F> map_impl<F,2> map2 (F const & f) { return map_impl<F,2>(f,true); }
template<class F> map_impl<F,1> map (F f) { return {std::move(f),true}; }
template<class F> map_impl<F,2> map2 (F f) { return {std::move(f),true}; }
// ----------- implementation ------------------------------------- // ----------- implementation -------------------------------------
template<typename F, int arity, bool is_vec, typename ... A> struct map_impl_result;
template<typename F, bool is_vec, typename A> struct map_impl_result<F,1,is_vec,A> {
typedef typename std::result_of<F(typename std::remove_reference<A>::type::value_type)>::type value_type;
typedef typename std::remove_reference<A>::type::domain_type domain_type;
F f;
typename std::add_const<A>::type a; // A is a T or a T& : add const to them.
domain_type domain() const { return a.domain(); }
template<typename ... Args> value_type operator() (Args && ... args) const { return f(a(std::forward<Args>(args)...)); }
friend std::ostream & operator<<(std::ostream & out, map_impl_result const & x){ return out<<"mapping result";}
// rest is only for vector
template<bool vec = is_vec>
TYPE_ENABLE_IFC(size_t,vec) size() const { return a.size();}
template<typename Args, bool vec=is_vec>
TYPE_ENABLE_IFC(value_type,vec) operator[] (Args && args) const { return f(a[std::forward<Args>(args)]);}
};
// possible to generalize to N order using tuple techniques ...
template<typename F, bool is_vec, typename A, typename B> struct map_impl_result<F,2,is_vec,A,B> {
typedef typename std::result_of<F(typename remove_cv_ref<A>::type::value_type, typename remove_cv_ref<B>::type::value_type)>::type value_type;
typedef typename remove_cv_ref<A>::type::domain_type domain_type;
F f;
typename std::add_const<A>::type a;
typename std::add_const<B>::type b;
domain_type domain() const { return a.domain(); }
template<typename ... Args> value_type operator() (Args && ... args) const { return f(a(std::forward<Args>(args)...),b(std::forward<Args>(args)...)); }
friend std::ostream & operator<<(std::ostream & out, map_impl_result const & x){ return out<<"mapping result";}
// rest is only for vector
template<bool vec = is_vec>
TYPE_ENABLE_IFC(size_t,vec) size() const { return a.size();}
template<typename Args, bool vec=is_vec>
TYPE_ENABLE_IFC(value_type,vec) operator[] (Args && args) const { return f(a[std::forward<Args>(args)],b[std::forward<Args>(args)]);}
};
template<typename ... T> struct _and;
template<typename T0, typename ... T> struct _and<T0, T...> : std::integral_constant<bool, T0::value && _and<T...>::value>{};
template<typename T> struct _and<T> : T{};
template<typename F, int arity, bool b, typename ... A> struct ImmutableCuboidArray<map_impl_result<F,arity,b,A...>> : std::true_type{};
template<typename F, int arity, bool b, typename ... A> struct ImmutableArray <map_impl_result<F,arity,b,A...>> : _and<typename ImmutableArray <typename std::remove_reference<A>::type>::type...>{};
template<typename F, int arity, bool b, typename ... A> struct ImmutableMatrix<map_impl_result<F,arity,b,A...>> : _and<typename ImmutableMatrix<typename std::remove_reference<A>::type>::type...>{};
template<typename F, int arity, bool b, typename ... A> struct ImmutableVector<map_impl_result<F,arity,b,A...>> : _and<typename ImmutableVector<typename std::remove_reference<A>::type>::type...>{};
//template<typename F, int arity, bool b, typename A> struct ImmutableArray <map_impl_result<F,arity,b,A>> : ImmutableArray <A>{};
//template<typename F, int arity, bool b, typename A> struct ImmutableMatrix<map_impl_result<F,arity,b,A>> : ImmutableMatrix<A>{};
//template<typename F, int arity, bool b, typename A> struct ImmutableVector<map_impl_result<F,arity,b,A>> : ImmutableVector<A>{};
// NB The bool is to make constructor not ambiguous // NB The bool is to make constructor not ambiguous
// clang on os X with lib++ has a pb otherwise (not clear what the pb is) // clang on os X with lib++ has a pb otherwise (not clear what the pb is)
template<class F> class map_impl<F,1> { template<class F, int arity> class map_impl {
F f; F f;
public : public :
// map_impl(F const & f_):f(f_) {} map_impl(F f_, bool):f(std::move(f_)) {}
map_impl(F const & f_, bool):f(f_) {}
map_impl(F && f_, bool):f(f_) {}
map_impl(map_impl const &) = default; map_impl(map_impl const &) = default;
//map_impl(map_impl &&) = default; map_impl(map_impl &&) = default;
map_impl & operator = (map_impl const &) = default;
map_impl & operator = (map_impl &&) = default;
template<typename A, typename Enable = void> class m_result; template<typename ... A> map_impl_result<F,arity,_and<typename ImmutableVector<A>::type...>::value,A...>
operator()(A&&... a) const { return {f,std::forward<A>(a)...}; }
template<typename A> class m_result<A,typename boost::enable_if<ImmutableArray<A> >::type> : TRIQS_CONCEPT_TAG_NAME(ImmutableArray) {
public:
typedef typename std::result_of<F(typename A::value_type)>::type value_type;
typedef typename A::domain_type domain_type;
A const & a; F f;
m_result(F const & f_, A const & a_):a(a_),f(f_) {}
domain_type domain() const { return a.domain(); }
template<typename ... Args> value_type operator() (Args const & ... args) const { return f(a(args...)); }
friend std::ostream & operator<<(std::ostream & out, m_result const & x){ return out<<"lazy matrix resulting of a mapping";}
};
template<typename A> class m_result<A,typename boost::enable_if<ImmutableMatrix<A> >::type> : TRIQS_CONCEPT_TAG_NAME(ImmutableMatrix) {
public:
typedef typename std::result_of<F(typename A::value_type)>::type value_type;
typedef typename A::domain_type domain_type;
A const & a; F f;
m_result(F const & f_, A const & a_):a(a_),f(f_) {}
domain_type domain() const { return a.domain(); }
template<typename ... Args> value_type operator() (Args const & ... args) const { return f(a(args...)); }
friend std::ostream & operator<<(std::ostream & out, m_result const & x){ return out<<"lazy matrix resulting of a mapping";}
};
template<typename A> class m_result<A,typename boost::enable_if<ImmutableVector<A> >::type> : TRIQS_CONCEPT_TAG_NAME(ImmutableVector) {
public:
typedef typename std::result_of<F(typename A::value_type)>::type value_type;
typedef typename A::domain_type domain_type;
A const & a; F f;
m_result(F const & f_, A const & a_):a(a_),f(f_) {}
domain_type domain() const { return a.domain(); }
size_t size() const { return a.size();}
template<typename ... Args> value_type operator() (Args const & ... args) const { return f(a(args...)); }
value_type operator() ( size_t i) const { return f(a(i)); }
friend std::ostream & operator<<(std::ostream & out, m_result const & x){ return out<<"lazy matrix resulting of a mapping";}
};
template<typename Sig> struct result;
template<typename This, typename A> struct result<This(A)> { typedef m_result<A> type;};
template< class A > m_result<A> operator()(A const & a) const {
//static_assert( (ImmutableCuboidArray<A>::value), "map : A does not model ImmutableCuboidArray");
return m_result<A>(f,a);
}
friend std::ostream & operator<<(std::ostream & out, map_impl const & x){ return out<<"map("<<"F"<<")";} friend std::ostream & operator<<(std::ostream & out, map_impl const & x){ return out<<"map("<<"F"<<")";}
}; };
// ----------- // TO DO : make this with preprocessor ....
template<class F> class map_impl<F,2> {
F f;
public :
map_impl(F const & f_, bool):f(f_) {}
map_impl(F && f_, bool):f(f_) {}
map_impl(map_impl const &) = default;
//map_impl(map_impl &&) = default;
template<class A, class B> class m_result : TRIQS_CONCEPT_TAG_NAME(ImmutableCuboidArray) {
static_assert( (std::is_same<typename A::domain_type, typename B::domain_type>::value), "type mismatch");
public:
typedef typename std::result_of<F(typename A::value_type,typename B::value_type)>::type value_type;
typedef typename A::domain_type domain_type;
A const & a; B const & b; F f;
m_result(F const & f_, A const & a_, B const & b_):a(a_),b(b_),f(f_) {
if (a.domain() != b.domain()) TRIQS_RUNTIME_ERROR<<"map2 : domain mismatch";
}
domain_type domain() const { return a.domain(); }
template<typename ... Args> value_type operator() (Args const & ... args) const { return f(a(args...),b(args...)); }
friend std::ostream & operator<<(std::ostream & out, m_result const & x){ return out<<"lazy matrix resulting of a mapping";}
};
template<typename Sig> struct result;
template<typename This, typename A,typename B> struct result<This(A,B)> { typedef m_result<A,B> type;};
template< class A, class B > m_result<A,B> operator()(A const & a, B const & b) {
static_assert( (ImmutableCuboidArray<A>::value), "map1 : A does not model ImmutableCuboidArray");
static_assert( (ImmutableCuboidArray<B>::value), "map1 : B does not model ImmutableCuboidArray");
return m_result<A,B>(f,a,b);
}
friend std::ostream & operator<<(std::ostream & out, map_impl const & x){ return out<<"map("<<"F"<<")";}
};
}}//namespace triqs::arrays }}//namespace triqs::arrays
#endif #endif

86
triqs/arrays/mapped_functions.hpp
View File

@ -22,84 +22,52 @@
#define TRIQS_ARRAYS_MAPPED_FNT_H #define TRIQS_ARRAYS_MAPPED_FNT_H
#include "./functional/map.hpp" #include "./functional/map.hpp"
#include <boost/preprocessor/seq/for_each.hpp> #include <boost/preprocessor/seq/for_each.hpp>
#include <boost/utility/result_of.hpp>
namespace triqs { namespace arrays { namespace triqs { namespace arrays {
using std::abs; //C++14 will simply be ...
struct abs_wrap { //template <typename A> decltype(auto) abs(A && a) { return map( [](auto const &x) { using std::abs; return abs(a);}, std::forward<A>(a));}
template<typename Sig> struct result;
template<typename This, typename A> struct result<This(A)> { typedef A type;};
template<typename This, typename A> struct result<This(std::complex<A>)> { typedef A type;};
static const int arity =1;
template<typename A> typename result<abs_wrap(A)>::type operator()( A const & a) const { return abs(a);}
};
template <typename A>
typename std::result_of<map_impl<abs_wrap>(A)>::type abs(A const & a) { return map(abs_wrap())(a); }
using std::pow; using std::pow;
template<typename T> TYPE_ENABLE_IF(T, std::is_integral<T>) pow(T x, int n) { return (n==0 ? 1 : pow(x,n-1)*x);} template<typename T> TYPE_ENABLE_IF(T, std::is_integral<T>) pow(T x, int n) { return (n==0 ? 1 : pow(x,n-1)*x);}
struct pow_wrap { struct pow_wrap {
template<typename Sig> struct result; int n;
template<typename This, typename A> struct result<This(A)> { typedef A type;}; template<typename A> A operator()( A const & a) const { return pow(a,n);}
static const int arity =1;
int n; pow_wrap(int n_): n(n_){}
template<typename A> typename result<abs_wrap(A)>::type operator()( A const & a) const { return pow(a,n);}
}; };
template <typename A> template <typename A>
typename boost::result_of<map_impl<pow_wrap>(A)>::type pow(A const & a, int n) { return map(pow_wrap(n))(a); } typename boost::lazy_enable_if_c<!ImmutableMatrix<A>::value,std::result_of<map_impl<pow_wrap,1>(A)>>::type
pow(A && a, int n) { return map(pow_wrap{n})(std::forward<A>(a)); }
#define MAP_IT(FNT) \ #define MAP_IT(FNT) \
using std::FNT;\ using std::FNT;\
struct FNT##_wrap {\ struct __triqs_##FNT##_wrap { template<typename A> auto operator()(A const & a) const DECL_AND_RETURN(FNT(a)); };\
template<typename Sig> struct result;\ template <typename A> \
template<typename This, typename A> struct result<This(A)> { typedef A type;};\ typename boost::lazy_enable_if_c<ImmutableCuboidArray<A>::value,std::result_of<map_impl<__triqs_##FNT##_wrap,1>(A)>>::type\
static const int arity =1;\ FNT(A && a) { return map(__triqs_##FNT##_wrap{})(std::forward<A>(a)); }
template<typename A> typename result<FNT##_wrap(A)>::type operator()( A const & a) const { return FNT(a);}\
};\ #define TRIQS_ARRAYS_MATH_FNT (abs)(real)(imag)(floor)
template <typename A>\
typename boost::result_of<map_impl<FNT##_wrap>(A)>::type FNT(A const & a) { return map(FNT##_wrap())(a); }
#define TRIQS_ARRAYS_MATH_FNT1 (cos)(sin)(tan)(cosh)(sinh)(tanh)(acos)(asin)(atan)(exp)(log)(sqrt)(floor)
#define AUX(r, data, elem) MAP_IT(elem) #define AUX(r, data, elem) MAP_IT(elem)
BOOST_PP_SEQ_FOR_EACH(AUX , nil , TRIQS_ARRAYS_MATH_FNT1); BOOST_PP_SEQ_FOR_EACH(AUX , nil , TRIQS_ARRAYS_MATH_FNT);
#undef AUX #undef AUX
#undef MAP_IT #undef MAP_IT
#undef TRIQS_ARRAYS_MATH_FNT
// Functions only defined for vector, array but NOT for matrix. They act element wise.
using std::real; #define MAP_IT(FNT) \
struct real_wrap { using std::FNT;\
template<typename Sig> struct result; struct __triqs_##FNT##_wrap { template<typename A> auto operator()(A const & a) const DECL_AND_RETURN(FNT(a)); };\
template<typename This, typename A> struct result<This(A)> { typedef A type;}; template <typename A> \
template<typename This, typename A> struct result<This(std::complex<A>)> { typedef A type;}; typename boost::lazy_enable_if_c<ImmutableArray<A>::value||ImmutableVector<A>::value,std::result_of<map_impl<__triqs_##FNT##_wrap,1>(A)>>::type\
static const int arity =1; FNT(A && a) { return map(__triqs_##FNT##_wrap{})(std::forward<A>(a)); }
// template<typename A> auto operator()( A const & a) const -> decltype(real(a)) { return real(a);}
template<typename A> typename result<real_wrap(A)>::type operator()( A const & a) const { return real(a);}
};
template <typename A>
//auto real(A const & a) -> decltype( map(real_wrap())(a)) { return map(real_wrap())(a); }
// typename std::result_of<map(real_wrap)(A)>::type real(A const & a) { return map(real_wrap())(a); }
typename boost::result_of<map_impl<real_wrap>(A)>::type real(A const & a) { return map(real_wrap())(a); }
using std::imag; #define TRIQS_ARRAYS_MATH_FNT (exp)(cos)(sin)(tan)(cosh)(sinh)(tanh)(acos)(asin)(atan)(log)(sqrt)
struct imag_wrap {
template<typename Sig> struct result;
template<typename This, typename A> struct result<This(A)> { typedef A type;};
template<typename This, typename A> struct result<This(std::complex<A>)> { typedef A type;};
static const int arity =1;
// template<typename A> auto operator()( A const & a) const -> decltype(imag(a)) { return imag(a);}
template<typename A> typename result<imag_wrap(A)>::type operator()( A const & a) const { return imag(a);}
};
template <typename A>
//auto imag(A const & a) -> decltype( map(imag_wrap())(a)) { return map(imag_wrap())(a); }
// typename std::result_of<map(imag_wrap)(A)>::type imag(A const & a) { return map(imag_wrap())(a); }
typename boost::result_of<map_impl<imag_wrap>(A)>::type imag(A const & a) { return map(imag_wrap())(a); }
#define AUX(r, data, elem) MAP_IT(elem)
BOOST_PP_SEQ_FOR_EACH(AUX , nil , TRIQS_ARRAYS_MATH_FNT);
#undef AUX
#undef MAP_IT
#undef TRIQS_ARRAYS_MATH_FNT
}}//namespace triqs::arrays }}//namespace triqs::arrays
#endif #endif

9
triqs/utility/concept_tools.hpp
View File

@ -31,13 +31,18 @@
#define TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT(MyBeautifulConcept) \ #define TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT(MyBeautifulConcept) \
struct TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept) {};\ struct TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept) {};\
template<typename T> struct MyBeautifulConcept : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept) , T> {}; template<typename T> struct MyBeautifulConcept : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept) , T> {};\
template<typename T> struct MyBeautifulConcept<T&> : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept), T> {};\
template<typename T> struct MyBeautifulConcept<T&&> : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept), T> {};
#define TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT_R_AUX(r, data, i, elem) BOOST_PP_COMMA_IF(i) TRIQS_CONCEPT_TAG_NAME(elem) #define TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT_R_AUX(r, data, i, elem) BOOST_PP_COMMA_IF(i) TRIQS_CONCEPT_TAG_NAME(elem)
#define TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT_R(MyBeautifulConcept,Rs) \ #define TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT_R(MyBeautifulConcept,Rs) \
struct TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept) : BOOST_PP_SEQ_FOR_EACH_I (TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT_R_AUX,nil,Rs) {};\ struct TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept) : BOOST_PP_SEQ_FOR_EACH_I (TRIQS_DEFINE_CONCEPT_AND_ASSOCIATED_TRAIT_R_AUX,nil,Rs) {};\
template<typename T> struct MyBeautifulConcept : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept), T> {}; template<typename T> struct MyBeautifulConcept : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept), T> {};\
template<typename T> struct MyBeautifulConcept<T&> : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept), T> {};\
template<typename T> struct MyBeautifulConcept<T&&> : boost::is_base_of<TRIQS_CONCEPT_TAG_NAME(MyBeautifulConcept), T> {};
#ifdef TRIQS_COMPILE_TIME_DEBUG #ifdef TRIQS_COMPILE_TIME_DEBUG
#define TRIQS_ASSERT_MODEL_CONCEPT(MyBeautifulConcept,T) BOOST_CONCEPT_ASSERT((BCC_##MyBeautifulConcept<T>)); #define TRIQS_ASSERT_MODEL_CONCEPT(MyBeautifulConcept,T) BOOST_CONCEPT_ASSERT((BCC_##MyBeautifulConcept<T>));