29 #if !defined( FEELPP_VF_VAL_HPP )
30 #define FEELPP_VF_VAL_HPP 1
32 # include <boost/preprocessor/comparison/less.hpp>
33 # include <boost/preprocessor/logical/and.hpp>
34 # include <boost/preprocessor/control/if.hpp>
35 # include <boost/preprocessor/list/at.hpp>
36 # include <boost/preprocessor/list/cat.hpp>
37 # include <boost/preprocessor/list/for_each_product.hpp>
38 # include <boost/preprocessor/logical/or.hpp>
39 # include <boost/preprocessor/tuple/to_list.hpp>
40 # include <boost/preprocessor/tuple/eat.hpp>
41 # include <boost/preprocessor/facilities/empty.hpp>
42 # include <boost/preprocessor/punctuation/comma.hpp>
43 # include <boost/preprocessor/facilities/identity.hpp>
46 #include <feel/feelcore/traits.hpp>
50 #if defined( FEELPP_HAS_QD_H ) && defined(FEELPP_HAS_MPFR)
51 # define VF_CHECK_ARITHMETIC_TYPE() \
52 BOOST_STATIC_ASSERT( (::boost::is_arithmetic<value_1_type>::value || \
53 ::boost::is_same<value_1_type, std::complex<float> >::value || \
54 ::boost::is_same<value_1_type, std::complex<double> >::value || \
55 ::boost::is_same<value_1_type,mp_type>::value || \
56 ::boost::is_same<value_1_type,dd_real>::value || \
57 ::boost::is_same<value_1_type,qd_real>::value) ); \
59 #elif defined( FEELPP_HAS_QD_H )
60 # define VF_CHECK_ARITHMETIC_TYPE() \
61 BOOST_STATIC_ASSERT( (::boost::is_arithmetic<value_1_type>::value || \
62 ::boost::is_same<value_1_type, std::complex<float> >::value || \
63 ::boost::is_same<value_1_type, std::complex<double> >::value || \
64 ::boost::is_same<value_1_type,dd_real>::value || \
65 ::boost::is_same<value_1_type,qd_real>::value) ); \
67 #elif defined( FEELPP_HAS_MPFR )
68 # define VF_CHECK_ARITHMETIC_TYPE() \
69 BOOST_STATIC_ASSERT( (::boost::is_arithmetic<value_1_type>::value || \
70 ::boost::is_same<value_1_type, std::complex<float> >::value || \
71 ::boost::is_same<value_1_type, std::complex<double> >::value || \
72 ::boost::is_same<value_1_type,mp_type>::value) ); \
75 # define VF_CHECK_ARITHMETIC_TYPE() \
76 BOOST_STATIC_ASSERT( ( ::boost::is_arithmetic<value_1_type>::value || \
77 ::boost::is_same<value_1_type, std::complex<float> >::value || \
78 ::boost::is_same<value_1_type, std::complex<double> >::value ) \
88 template <
typename ExprT1 >
91 public UnaryFunctor<typename ExprT1::value_type>
95 static const size_type context = ExprT1::context;
96 static const bool is_terminal = ExprT1::is_terminal;
98 static const uint16_type imorder = ExprT1::imorder;
99 static const bool imIsPoly = ExprT1::imIsPoly;
101 template<
typename Func>
102 struct HasTestFunction
104 static const bool result =
false;
107 template<
typename Func>
108 struct HasTrialFunction
110 static const bool result =
false;
113 typedef UnaryFunctor<typename ExprT1::value_type> super;
114 typedef typename super::functordomain_type functordomain_type;
115 typedef typename super::functordomain_ptrtype functordomain_ptrtype;
116 typedef ExprT1 expression_1_type;
117 typedef Val<ExprT1> this_type;
118 typedef typename expression_1_type::value_type value_1_type;
119 typedef value_1_type value_type;
121 VF_CHECK_ARITHMETIC_TYPE()
123 explicit Val( expression_1_type const& __expr1 )
125 super( "value", functordomain_ptrtype( new UnboundedDomain<value_type>() ) ),
128 DVLOG(2) <<
"Val::Val default constructorn";
131 Val( Val
const& __vfp )
133 super(
"value", functordomain_ptrtype( new UnboundedDomain<value_type>() ) ),
134 M_expr_1( __vfp.M_expr_1 )
136 DVLOG(2) <<
"Val::Val copy constructorn";
139 bool isSymetric()
const
144 void eval(
int nx, value_type
const* x, value_type* f )
const
146 for (
int i = 0; i < nx; ++i )
150 expression_1_type
const& expression()
const
155 template<
typename Geo_t,
typename Basis_i_t,
typename Basis_j_t = Basis_i_t>
158 typedef this_type expression_type;
160 typedef typename expression_1_type::template tensor<Geo_t> tensor2_expr_type;
161 typedef typename tensor2_expr_type::value_type value_type;
162 typedef typename mpl::if_<fusion::result_of::has_key<Geo_t,vf::detail::gmc<0> >,
163 mpl::identity<vf::detail::gmc<0> >,
164 mpl::identity<vf::detail::gmc<1> > >::type::type key_type;
165 typedef typename fusion::result_of::value_at_key<Geo_t,key_type>::type::element_type* gmc_ptrtype;
166 typedef typename fusion::result_of::value_at_key<Geo_t,key_type>::type::element_type gmc_type;
167 typedef typename tensor2_expr_type::shape shape;
171 static const bool value = tensor2_expr_type::is_zero::value;
174 template<
typename ExprT>
175 tensor( ExprT
const& expr, Geo_t
const& geom, Basis_i_t
const& , Basis_j_t
const& )
177 M_expr( expr.expression(), geom ),
178 M_gmc( fusion::at_key<key_type>( geom ).get() ),
179 M_loc( boost::extents[M_gmc->nPoints()][shape::M][shape::N] )
183 template<
typename ExprT>
184 tensor( ExprT
const& expr,Geo_t
const& geom, Basis_i_t
const& )
186 M_expr( expr.expression(), geom ),
187 M_gmc( fusion::at_key<key_type>( geom ).get() ),
188 M_loc( boost::extents[M_gmc->nPoints()][shape::M][shape::N] )
192 template<
typename ExprT>
193 tensor( ExprT
const& expr, Geo_t
const& geom )
195 M_expr( expr.expression(), geom ),
196 M_gmc( fusion::at_key<key_type>( geom ).get() ),
197 M_loc( boost::extents[M_gmc->nPoints()][shape::M][shape::N] )
201 template<
typename IM>
202 void init( IM
const& im )
206 void update( Geo_t
const& geom, Basis_i_t
const& , Basis_j_t
const& )
210 void update( Geo_t
const& geom, Basis_i_t
const& )
214 void update( Geo_t
const& geom )
216 M_expr.update( geom );
218 for (
int q = 0; q < M_gmc->nPoints(); ++q )
219 for (
int c1 = 0; c1 < shape::M; ++c1 )
220 for (
int c2 = 0; c2 < shape::N; ++c2 )
222 M_loc[q][c1][c2] = M_expr.evalq( c1, c2, q );
225 void update( Geo_t
const& geom, uint16_type face )
227 M_expr.update( geom, face );
229 for (
int q = 0; q < M_gmc->nPoints(); ++q )
230 for (
int c1 = 0; c1 < shape::M; ++c1 )
231 for (
int c2 = 0; c2 < shape::N; ++c2 )
233 M_loc[q][c1][c2] = M_expr.evalq( c1, c2, q );
238 evalijq( uint16_type , uint16_type , uint16_type c1, uint16_type c2, uint16_type q )
const
240 return evalq( c1, c2, q );
242 template<
int PatternContext>
244 evalijq( uint16_type , uint16_type , uint16_type c1, uint16_type c2, uint16_type q,
245 mpl::int_<PatternContext> )
const
247 return evalq( c1, c2, q );
251 evaliq( uint16_type , uint16_type c1, uint16_type c2, uint16_type q )
const
253 return evalq( c1, c2, q );
256 evalq( uint16_type c1, uint16_type c2, uint16_type q )
const
258 return evalq( c1, c2, q, mpl::int_<shape::rank>() );
262 evalq( uint16_type c1, uint16_type c2, uint16_type q, mpl::int_<0> )
const
264 Feel::detail::ignore_unused_variable_warning( c1 );
265 Feel::detail::ignore_unused_variable_warning( c2 );
266 return M_loc[q][0][0];
269 evalq( uint16_type c1, uint16_type c2, uint16_type q, mpl::int_<1> )
const
271 if ( shape::M > shape::N )
272 return M_loc[q][c1][0];
274 return M_loc[q][0][c2];
277 evalq( uint16_type c1, uint16_type c2, uint16_type q, mpl::int_<2> )
const
279 return M_loc[q][c1][c2];
282 tensor2_expr_type M_expr;
284 boost::multi_array<value_type,3> M_loc;
290 expression_1_type M_expr_1;
299 template<
typename ExprT1>
301 Expr< Val<typename mpl::if_<boost::is_arithmetic<ExprT1>,
302 mpl::identity<Cst<ExprT1> >,
303 mpl::identity<ExprT1> >::type::type > >
306 typedef typename mpl::if_<boost::is_arithmetic<ExprT1>,
307 mpl::identity<Cst<ExprT1> >,
308 mpl::identity<ExprT1> >::type::type t1;
309 typedef Val<t1> expr_t;
310 return Expr< expr_t >( expr_t( t1( __e1 ) ) );