From 4146827dd6622db6d5051a2cfaf6cd3cfe23d146 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 15 Sep 2010 02:44:35 +0000 Subject: [PATCH] Make compile with gcc 3.3.x. It has trouble finding the proper overload of std::sqrt in some cases. git-svn-id: https://svn.dealii.org/trunk@21971 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/source/polynomial.cc | 108 +++++++++++++++--------------- 1 file changed, 54 insertions(+), 54 deletions(-) diff --git a/deal.II/base/source/polynomial.cc b/deal.II/base/source/polynomial.cc index e82e048956..a84e6379b3 100644 --- a/deal.II/base/source/polynomial.cc +++ b/deal.II/base/source/polynomial.cc @@ -1,8 +1,8 @@ //--------------------------------------------------------------------------- -// $Id$ +// $Id$ // Version: $Name$ // -// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2010 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -31,7 +31,7 @@ DEAL_II_NAMESPACE_OPEN // to be a problem since we only need it on very rare occasions. if // someone finds this is a bottleneck, feel free to replace it by a // more fine-grained solution -namespace +namespace { Threads::ThreadMutex coefficients_lock; } @@ -60,8 +60,8 @@ namespace Polynomials Assert (coefficients.size() > 0, ExcEmptyObject()); Assert (values.size() > 0, ExcZero()); const unsigned int values_size=values.size(); - - + + // if we only need the value, then // call the other function since // that is significantly faster @@ -116,7 +116,7 @@ namespace Polynomials { *c *= f; f *= factor; - } + } } @@ -152,7 +152,7 @@ namespace Polynomials return *this; } - + template Polynomial& Polynomial::operator *= (const Polynomial& p) @@ -161,16 +161,16 @@ namespace Polynomials unsigned int new_degree = this->degree() + p.degree(); std::vector new_coefficients(new_degree+1, 0.); - + for (unsigned int i=0; icoefficients.size(); ++j) new_coefficients[i+j] += this->coefficients[j]*p.coefficients[i]; this->coefficients = new_coefficients; - + return *this; } - + template Polynomial& Polynomial::operator += (const Polynomial& p) @@ -186,7 +186,7 @@ namespace Polynomials return *this; } - + template Polynomial& Polynomial::operator -= (const Polynomial& p) @@ -202,13 +202,13 @@ namespace Polynomials return *this; } - + template template void Polynomial::shift(std::vector& coefficients, const number2 offset) - { + { #ifdef DEAL_II_LONG_DOUBLE_LOOP_BUG AssertThrow (false, ExcMessage("Sorry, but the compiler you are using has a bug that disallows " @@ -220,12 +220,12 @@ namespace Polynomials // args. note that this code is // actually unreachable coefficients[0] = offset; -#else +#else // Copy coefficients to a vector of // accuracy given by the argument std::vector new_coefficients(coefficients.begin(), coefficients.end()); - + // Traverse all coefficients from // c_1. c_0 will be modified by // higher degrees, only. @@ -242,7 +242,7 @@ namespace Polynomials // needed and computed // successively. number2 offset_power = offset; - + // Compute (x+offset)^d // and modify all values c_k // with k Polynomial Polynomial::derivative () const @@ -297,7 +297,7 @@ namespace Polynomials return Polynomial (newcoefficients); } - + template Polynomial @@ -310,7 +310,7 @@ namespace Polynomials return Polynomial (newcoefficients); } - + template void @@ -335,15 +335,15 @@ namespace Polynomials result[n] = coefficient; return result; } - - + + template Monomial::Monomial (unsigned int n, double coefficient) : Polynomial(make_vector(n, coefficient)) {} - - + + template std::vector > Monomial::generate_complete_basis (const unsigned int degree) @@ -364,7 +364,7 @@ namespace Polynomials - std::vector + std::vector LagrangeEquidistant::compute_coefficients (const unsigned int n, const unsigned int support_point) { @@ -375,7 +375,7 @@ namespace Polynomials Assert(support_point > Lagrange::generate_complete_basis (const std::vector >& points) { @@ -621,7 +621,7 @@ namespace Polynomials std::vector linear(2, 1.); // We start with a constant polynomial std::vector one(1, 1.); - + for (unsigned int i=0;i(*c0); c1 = new std::vector(*c1); - + Polynomial::shift (*c0, -1.); Polynomial::scale(*c0, 2.); Polynomial::shift (*c1, -1.); @@ -776,11 +776,11 @@ namespace Polynomials coefficients_lock.acquire (); std::vector *ck = new std::vector(k+1); - + const double a = 1./(k); const double b = a*(2*k-1); const double c = a*(k-1); - + (*ck)[k] = b*(*recursive_coefficients[k-1])[k-1]; (*ck)[k-1] = b*(*recursive_coefficients[k-1])[k-2]; for (unsigned int i=1 ; i<= k-2 ; ++i) @@ -862,8 +862,8 @@ std::vector Lobatto::compute_coefficients (const unsigned int p) { std::vector coefficients (3); coefficients[0] = 0.0; - coefficients[1] = -1.0 * std::sqrt (3); - coefficients[2] = std::sqrt (3); + coefficients[1] = -1.0 * std::sqrt (3.); + coefficients[2] = std::sqrt (3.); return coefficients; } @@ -872,8 +872,8 @@ std::vector Lobatto::compute_coefficients (const unsigned int p) { std::vector legendre_coefficients_tmp1 (p); std::vector legendre_coefficients_tmp2 (p - 1); - coefficients[0] = -1.0 * std::sqrt (3); - coefficients[1] = 2.0 * std::sqrt (3); + coefficients[0] = -1.0 * std::sqrt (3.); + coefficients[1] = 2.0 * std::sqrt (3.); legendre_coefficients_tmp1[0] = 1.0; for (unsigned int i = 2; i < p; ++i) { @@ -883,13 +883,13 @@ std::vector Lobatto::compute_coefficients (const unsigned int p) { for (unsigned int j = 0; j < i; ++j) legendre_coefficients_tmp1[j] = coefficients[j]; - coefficients[0] = std::sqrt (2 * i + 1) * ((1.0 - 2 * i) * legendre_coefficients_tmp1[0] / std::sqrt (2 * i - 1) + (1.0 - i) * legendre_coefficients_tmp2[0] / std::sqrt (2 * i - 3)) / i; + coefficients[0] = std::sqrt (2 * i + 1.) * ((1.0 - 2 * i) * legendre_coefficients_tmp1[0] / std::sqrt (2 * i - 1.) + (1.0 - i) * legendre_coefficients_tmp2[0] / std::sqrt (2 * i - 3.)) / i; for (unsigned int j = 1; j < i - 1; ++j) - coefficients[j] = std::sqrt (2 * i + 1) * (std::sqrt (2 * i - 1) * (2.0 * legendre_coefficients_tmp1[j - 1] - legendre_coefficients_tmp1[j]) + (1.0 - i) * legendre_coefficients_tmp2[j] / std::sqrt (2 * i - 3)) / i; + coefficients[j] = std::sqrt (2 * i + 1.) * (std::sqrt (2 * i - 1.) * (2.0 * legendre_coefficients_tmp1[j - 1] - legendre_coefficients_tmp1[j]) + (1.0 - i) * legendre_coefficients_tmp2[j] / std::sqrt (2 * i - 3.)) / i; - coefficients[i - 1] = std::sqrt (4 * i * i - 1) * (2.0 * legendre_coefficients_tmp1[i - 2] - legendre_coefficients_tmp1[i - 1]) / i; - coefficients[i] = 2.0 * std::sqrt (4 * i * i - 1) * legendre_coefficients_tmp1[i - 1] / i; + coefficients[i - 1] = std::sqrt (4 * i * i - 1.) * (2.0 * legendre_coefficients_tmp1[i - 2] - legendre_coefficients_tmp1[i - 1]) / i; + coefficients[i] = 2.0 * std::sqrt (4 * i * i - 1.) * legendre_coefficients_tmp1[i - 1] / i; } for (int i = p; i > 0; --i) @@ -941,20 +941,20 @@ std::vector > Lobatto::generate_complete_basis (const unsigne // until we quit this function Threads::ThreadMutex::ScopedLock lock(coefficients_lock); - // The first 2 coefficients + // The first 2 coefficients // are hard-coded if (k==0) k=1; // check: does the information // already exist? if ( (recursive_coefficients.size() < k+1) || - ((recursive_coefficients.size() >= k+1) && + ((recursive_coefficients.size() >= k+1) && (recursive_coefficients[k] == 0)) ) // no, then generate the // respective coefficients { recursive_coefficients.resize (k+1, 0); - + if (k<=1) { // create coefficients @@ -991,7 +991,7 @@ std::vector > Lobatto::generate_complete_basis (const unsigne (*c2)[0] = 0.*a; (*c2)[1] = -4.*a; (*c2)[2] = 4.*a; - + recursive_coefficients[2] = c2; } else @@ -1009,15 +1009,15 @@ std::vector > Lobatto::generate_complete_basis (const unsigne coefficients_lock.acquire (); std::vector *ck = new std::vector(k+1); - + const double a = 1.; //1./(2.*k); (*ck)[0] = - a*(*recursive_coefficients[k-1])[0]; - + for (unsigned int i=1; i<=k-1; ++i) (*ck)[i] = a*( 2.*(*recursive_coefficients[k-1])[i-1] - (*recursive_coefficients[k-1])[i] ); - + (*ck)[k] = a*2.*(*recursive_coefficients[k-1])[k-1]; // for even degrees, we need // to add a multiple of @@ -1030,7 +1030,7 @@ std::vector > Lobatto::generate_complete_basis (const unsigne (*ck)[1] += b*(*recursive_coefficients[2])[1]; (*ck)[2] += b*(*recursive_coefficients[2])[2]; - } + } // finally assign the newly // created vector to the // const pointer in the -- 2.39.5