From 54407c0f47bd0734942456be2d18add9d5146d0f Mon Sep 17 00:00:00 2001 From: kronbichler Date: Tue, 2 Aug 2011 12:19:58 +0000 Subject: [PATCH] Fix scalar product with complex number: return real_type, not Number. Make unroll() work also for vectors that are not based on doubles by moving it to the .h file. This obviates tensor.cc git-svn-id: https://svn.dealii.org/trunk@23991 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/base/tensor.h | 85 +++++++++++++++++-- deal.II/include/deal.II/base/tensor_base.h | 98 ++++++++++++++++++---- deal.II/source/base/tensor.cc | 93 -------------------- 3 files changed, 159 insertions(+), 117 deletions(-) delete mode 100644 deal.II/source/base/tensor.cc diff --git a/deal.II/include/deal.II/base/tensor.h b/deal.II/include/deal.II/base/tensor.h index 8cfb1f71e6..effb6875a1 100644 --- a/deal.II/include/deal.II/base/tensor.h +++ b/deal.II/include/deal.II/base/tensor.h @@ -76,6 +76,22 @@ class Tensor */ typedef Tensor value_type; + /** + * Declare a type that has holds + * real-valued numbers with the same + * precision as the template argument to + * this class. For std::complex, + * this corresponds to type number, and + * it is equal to Number for all other + * cases. See also the respective field + * in Vector. + * + * This typedef is used to + * represent the return type of + * norms. + */ + typedef typename numbers::NumberTraits::real_type real_type; + /** * Declare an array type which * can be used to initialize an @@ -185,7 +201,7 @@ class Tensor * i.e. the square root of the sum of * squares of all entries. */ - Number norm () const; + real_type norm () const; /** * Return the square of the @@ -199,7 +215,7 @@ class Tensor * may also be useful in other * contexts. */ - Number norm_square () const; + real_type norm_square () const; /** * Fill a vector with all tensor elements. @@ -210,7 +226,8 @@ class Tensor * usual in C++, the rightmost * index of the tensor marches fastest. */ - void unroll (Vector & result) const; + template + void unroll (Vector & result) const; /** @@ -264,8 +281,9 @@ class Tensor /** * Help function for unroll. */ - void unroll_recursion(Vector &result, - unsigned int &start_index) const; + template + void unroll_recursion(Vector &result, + unsigned int &start_index) const; // make the following class a // friend to this class. in principle, @@ -295,6 +313,7 @@ Tensor::Tensor () } + template inline Tensor::Tensor (const array_type &initializer) @@ -304,6 +323,7 @@ Tensor::Tensor (const array_type &initializer) } + template inline typename Tensor::value_type& @@ -315,6 +335,7 @@ Tensor::operator[] (const unsigned int i) } + template inline const typename Tensor::value_type& @@ -326,6 +347,7 @@ Tensor::operator[] (const unsigned int i) const } + template inline Tensor & @@ -337,12 +359,13 @@ Tensor::operator = (const Tensor &t) } + template inline Tensor & Tensor::operator = (const Number d) { - Assert (d==0, ExcMessage ("Only assignment with zero is allowed")); + Assert (d==Number(0), ExcMessage ("Only assignment with zero is allowed")); for (unsigned int i=0; i::operator = (const Number d) } + template inline bool @@ -361,6 +385,7 @@ Tensor::operator == (const Tensor &p) const } + template inline bool @@ -370,6 +395,7 @@ Tensor::operator != (const Tensor &p) const } + template inline Tensor & @@ -381,6 +407,7 @@ Tensor::operator += (const Tensor &p) } + template inline Tensor & @@ -392,6 +419,7 @@ Tensor::operator -= (const Tensor &p) } + template inline Tensor & @@ -403,6 +431,7 @@ Tensor::operator *= (const Number s) } + template inline Tensor & @@ -414,6 +443,7 @@ Tensor::operator /= (const Number s) } + template inline Tensor @@ -428,6 +458,7 @@ Tensor::operator + (const Tensor &t) const } + template inline Tensor @@ -442,6 +473,7 @@ Tensor::operator - (const Tensor &t) const } + template inline Tensor @@ -456,19 +488,23 @@ Tensor::operator - () const } + template inline -Number Tensor::norm () const +typename Tensor::real_type +Tensor::norm () const { return std::sqrt (norm_square()); } + template inline -Number Tensor::norm_square () const +typename Tensor::real_type +Tensor::norm_square () const { - Number s = 0; + real_type s = 0; for (unsigned int i=0; i::norm_square () const } + +template +template +inline +void +Tensor::unroll (Vector &result) const +{ + AssertDimension (result.size(),std::pow(static_cast(dim),rank_)); + unsigned index = 0; + unroll_recursion (result, index); +} + + + +template +template +inline +void +Tensor::unroll_recursion (Vector &result, + unsigned int &index) const +{ + for (unsigned i=0; i inline void Tensor::clear () @@ -485,6 +550,7 @@ void Tensor::clear () } + template inline std::size_t @@ -494,6 +560,7 @@ Tensor::memory_consumption () } + template template inline diff --git a/deal.II/include/deal.II/base/tensor_base.h b/deal.II/include/deal.II/base/tensor_base.h index bbfc19eb9a..373be4f141 100644 --- a/deal.II/include/deal.II/base/tensor_base.h +++ b/deal.II/include/deal.II/base/tensor_base.h @@ -102,6 +102,22 @@ class Tensor<0,dim,Number> typedef Number value_type; + /** + * Declare a type that has holds + * real-valued numbers with the same + * precision as the template argument to + * this class. For std::complex, + * this corresponds to type number, and + * it is equal to Number for all other + * cases. See also the respective field + * in Vector. + * + * This typedef is used to + * represent the return type of + * norms. + */ + typedef typename numbers::NumberTraits::real_type real_type; + /** * Constructor. Set to zero. */ @@ -221,7 +237,7 @@ class Tensor<0,dim,Number> * l2 norm of * the vector. */ - Number norm () const; + real_type norm () const; /** * Return the square of the @@ -235,7 +251,7 @@ class Tensor<0,dim,Number> * may also be useful in other * contexts. */ - Number norm_square () const; + real_type norm_square () const; /** * Reset all values to zero. @@ -345,6 +361,22 @@ class Tensor<1,dim,Number> typedef Number value_type; + /** + * Declare a type that has holds + * real-valued numbers with the same + * precision as the template argument to + * this class. For std::complex, + * this corresponds to type number, and + * it is equal to Number for all other + * cases. See also the respective field + * in Vector. + * + * This typedef is used to + * represent the return type of + * norms. + */ + typedef typename numbers::NumberTraits::real_type real_type; + /** * Declare an array type which can * be used to initialize statically @@ -494,7 +526,7 @@ class Tensor<1,dim,Number> * l2 norm of * the vector. */ - Number norm () const; + real_type norm () const; /** * Return the square of the @@ -508,7 +540,7 @@ class Tensor<1,dim,Number> * may also be useful in other * contexts. */ - Number norm_square () const; + real_type norm_square () const; /** * Reset all values to zero. @@ -539,7 +571,8 @@ class Tensor<1,dim,Number> * usual in C++, the rightmost * index marches fastest. */ - void unroll (Vector &result) const; + template + void unroll (Vector &result) const; /** * Determine an estimate for @@ -596,8 +629,9 @@ class Tensor<1,dim,Number> * case, even if it should be * public for your compiler. */ - void unroll_recursion (Vector &result, - unsigned int &start_index) const; + template + void unroll_recursion (Vector &result, + unsigned int &start_index) const; private: /** @@ -810,18 +844,20 @@ Tensor<0,dim,Number> Tensor<0,dim,Number>::operator - () const template inline -Number Tensor<0,dim,Number>::norm () const +typename Tensor<0,dim,Number>::real_type +Tensor<0,dim,Number>::norm () const { - return std::abs (value); + return numbers::NumberTraits::abs (value); } template inline -Number Tensor<0,dim,Number>::norm_square () const +typename Tensor<0,dim,Number>::real_type +Tensor<0,dim,Number>::norm_square () const { - return value*value; + return numbers::NumberTraits::abs_square (value); } @@ -833,6 +869,8 @@ void Tensor<0,dim,Number>::clear () value = 0; } + + template template inline @@ -1125,7 +1163,8 @@ Tensor<1,dim,Number> Tensor<1,dim,Number>::operator - () const template inline -Number Tensor<1,dim,Number>::norm () const +typename Tensor<1,dim,Number>::real_type +Tensor<1,dim,Number>::norm () const { return std::sqrt (norm_square()); } @@ -1134,17 +1173,46 @@ Number Tensor<1,dim,Number>::norm () const template inline -Number Tensor<1,dim,Number>::norm_square () const +typename Tensor<1,dim,Number>::real_type +Tensor<1,dim,Number>::norm_square () const { - Number s = 0; + real_type s = 0; for (unsigned int i=0; i::abs_square(values[i]); return s; } +template +template +inline +void +Tensor<1,dim,Number>::unroll (Vector &result) const +{ + Assert (result.size()==dim, + ExcDimensionMismatch(dim, result.size())); + + unsigned index = 0; + unroll_recursion (result,index); +} + + + +template +template +inline +void +Tensor<1,dim,Number>::unroll_recursion (Vector &result, + unsigned int &index) const +{ + for (unsigned i=0; i inline void Tensor<1,dim,Number>::clear () diff --git a/deal.II/source/base/tensor.cc b/deal.II/source/base/tensor.cc deleted file mode 100644 index 5ade4a8554..0000000000 --- a/deal.II/source/base/tensor.cc +++ /dev/null @@ -1,93 +0,0 @@ -//--------------------------------------------------------------------------- -// $Id$ -// Version: $Name$ -// -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2005, 2006 by the deal.II authors -// -// This file is subject to QPL and may not be distributed -// without copyright and license information. Please refer -// to the file deal.II/doc/license.html for the text and -// further information on this license. -// -//--------------------------------------------------------------------------- - - -#include -#include -#include - -DEAL_II_NAMESPACE_OPEN - - -// storage for static variables -template -const unsigned int Tensor<1,dim,Number>::dimension; - -template -const unsigned int Tensor::dimension; - - -template -void -Tensor<1,dim,Number>::unroll (Vector &result) const -{ - Assert (result.size()==dim, - ExcDimensionMismatch(dim, result.size())); - - unsigned index = 0; - unroll_recursion (result,index); -} - - - -template -void -Tensor::unroll (Vector &result) const -{ - Assert(result.size()==std::pow(static_cast(dim),rank_), - ExcDimensionMismatch(static_cast(std::pow(static_cast(dim),rank_)), - result.size())); - - unsigned index = 0; - unroll_recursion (result, index); -} - - - -template -void -Tensor::unroll_recursion (Vector &result, - unsigned int &index) const -{ - for (unsigned i=0; i -void -Tensor<1,dim,Number>::unroll_recursion (Vector &result, - unsigned int &index) const -{ - for (unsigned i=0; i; -template class Tensor<1, 2, double>; -template class Tensor<1, 3, double>; -template class Tensor<2, 1, double>; -template class Tensor<2, 2, double>; -template class Tensor<2, 3, double>; -template class Tensor<3, 1, double>; -template class Tensor<3, 2, double>; -template class Tensor<3, 3, double>; -template class Tensor<4, 1, double>; -template class Tensor<4, 2, double>; -template class Tensor<4, 3, double>; - -DEAL_II_NAMESPACE_CLOSE -- 2.39.5