1 #ifndef DUNE_LOCALFUNCTIONS_UTILITY_FIELD_HH
2 #define DUNE_LOCALFUNCTIONS_UTILITY_FIELD_HH
5 #include <alglib/amp.h>
6 #warning ALGLIB support is deprecated and will be dropped after DUNE 2.2 (cf. FS#931)
9 #include <dune/common/gmpfield.hh>
10 #include <dune/common/fvector.hh>
11 #include <dune/common/fmatrix.hh>
29 template<
class Field >
32 operator Field ()
const
38 template<
class Field >
44 template<
class Field >
50 template<
class Field >
56 template<
class Field >
78 template<
class Field >
81 operator Field ()
const
91 template<
unsigned int precision >
92 struct Zero< amp::ampf< precision > >
94 typedef amp::ampf< precision > Field;
95 operator Field ()
const
106 template<
unsigned int precision >
107 struct Zero< GMPField< precision > >
109 typedef GMPField< precision > Field;
110 operator Field ()
const
121 template<
class Field >
127 template<
class Field >
133 template<
class Field >
134 inline bool operator< ( const Zero< Field > &,
const Field &f )
139 template<
class Field >
140 inline bool operator< ( const Field &f, const Zero< Field > & )
145 template<
class Field >
151 template<
class Field >
173 template<
class F2,
class F1 >
180 template<
unsigned int precision >
181 inline void field_cast (
const amp::ampf< precision > &f1,
double &f2 )
186 template<
unsigned int precision >
187 inline void field_cast (
const amp::ampf< precision > &f1,
long double &f2 )
194 template<
unsigned int precision >
195 inline void field_cast (
const Dune::GMPField< precision > &f1,
double &f2 )
200 template<
unsigned int precision >
201 inline void field_cast (
const Dune::GMPField< precision > &f1,
long double &f2 )
207 template<
class F2,
class F1,
int dim >
208 inline void field_cast (
const Dune::FieldVector< F1, dim > &f1, Dune::FieldVector< F2, dim > &f2 )
210 for(
int d = 0; d < dim; ++d )
213 template<
class F2,
class F1 >
214 inline void field_cast (
const Dune::FieldVector< F1, 1 > &f1, F2 &f2 )
218 template<
class F2,
class F1 >
219 inline void field_cast (
const F1 &f1, Dune::FieldVector< F2, 1 > &f2 )
224 template<
class F2,
class F1,
int rdim,
int cdim >
225 inline void field_cast (
const Dune::FieldMatrix< F1, rdim, cdim > &f1, Dune::FieldMatrix< F2, rdim, cdim > &f2 )
227 for(
int r = 0; r < rdim; ++r )
230 template<
class F2,
class F1 >
231 inline void field_cast (
const Dune::FieldMatrix<F1,1,1> &f1, Dune::FieldMatrix< F2, 1,1 > &f2 )
235 template<
class F2,
class F1 >
236 inline void field_cast (
const Dune::FieldMatrix< F1, 1,1 > &f1, F2 &f2 )
240 template<
class F2,
class F1 >
241 inline void field_cast (
const F1 &f1, Dune::FieldMatrix< F2, 1,1 > &f2 )
245 template<
class F2,
class F1 >
246 inline void field_cast (
const Dune::FieldVector<F1,1> &f1, Dune::FieldMatrix< F2, 1,1 > &f2 )
250 template<
class F2,
class F1 >
251 inline void field_cast (
const Dune::FieldMatrix<F1,1,1> &f1, Dune::FieldVector< F2, 1 > &f2 )
256 template<
class F2,
class F1 >
257 inline void field_cast (
const Dune::FieldVector< F1, 1 > &f1, Dune::FieldVector<F2, 1> &f2 )
262 template<
class F2,
class V >
267 template<
class F2,
class F1,
int dim >
270 typedef Dune::FieldVector<F2,dim>
type;
272 template<
class F2,
class F1,
int dim1,
int dim2>
275 typedef Dune::FieldMatrix<F2,dim1,dim2>
type;
277 template<
class F2,
class V >
292 template <
class Field>
298 static const unsigned int value = 64;
304 static const unsigned int value = 80;
310 static const unsigned int value = 32;
314 template<
unsigned int precision >
315 struct Precision< amp::ampf< precision > >
317 static const unsigned int value = precision;
321 template<
unsigned int precision >
322 struct Precision< GMPField< precision > >
324 static const unsigned int value = precision;
331 template <
class Field,
unsigned int sum>
338 template<
unsigned int precision,
unsigned int sum >
341 typedef amp::ampf<precision+sum>
Type;
345 template<
unsigned int precision,
unsigned int sum >
348 typedef GMPField<precision+sum>
Type;
358 template<
unsigned int precision >
361 const amp::ampf< precision > &
value )
363 return out << value.toDec();
366 template<
unsigned int precision >
367 inline amp::ampf< precision > sqrt (
const amp::ampf< precision > &a )
369 return amp::sqrt( a );
372 template<
unsigned int precision >
373 inline amp::ampf< precision > abs (
const amp::ampf< precision > &a )
375 return amp::abs( a );
377 #endif // #if HAVE_ALGLIB
384 namespace GenericGeometry
390 template<
class Field >
394 template<
unsigned int precision >
397 typedef amp::ampf< precision > Field;
399 static Field abs (
const Field &x ) {
return amp::abs( x ); }
401 #endif // #if HAVE_ALGLIB
407 #endif // #ifndef DUNE_LOCALFUNCTIONS_UTILITY_FIELD_HH