Fix for GCC 4.7
Kyrre Ness Sjøbæk <[email protected]> Mon, 27 Aug 2012 16:53:28 +0200
| Newsgroups | gmane.comp.lib.mtl.devel |
|---|---|
| Message-ID | <[email protected]> |
This is a multi-part message in MIME format.
--------------090106020405040407040902
Content-Type: text/plain; charset=ISO-8859-1; format=flowed
Content-Transfer-Encoding: 7bit
Hello,
I recently needed to compile a large simulation package containing ITL using GCC 4.7, and hit a few problems due to the new compiler version being stricter. To fix this, I changed a few lines, and it now compiles. The modified files are attached, and a diff is shown below.
Does this fix look correct to you? If so, please accept the patch.
Regards,
Kyrre Sjobak
> diff --git a/itl/itl/itl.h b/itl/itl/itl.h
> index 17fdb69..b91c68f 100644
> --- a/itl/itl/itl.h
> +++ b/itl/itl/itl.h
> @@ -226,7 +226,7 @@ class noisy_iteration : public basic_iteration<Real> {
>
> Real normr_ = std::abs(itl::two_norm(r));
> bool ret;
> - if (converged(normr_)) {
> + if (this->converged(normr_)) {
> ret = true;
> } else if (Super::i < Super::max_iter) {
> ret = false;
> @@ -244,7 +244,7 @@ class noisy_iteration : public basic_iteration<Real> {
> using std::endl;
>
> bool ret;
> - if (converged(r)) {
> + if (this->converged(r)) {
> ret = true;
> } else if (Super::i < Super::max_iter) {
> ret = false;
> @@ -263,7 +263,7 @@ class noisy_iteration : public basic_iteration<Real> {
> using std::endl;
>
> bool ret;
> - if (converged(std::abs(r))) {
> + if (this->converged(std::abs(r))) {
> ret = true;
> } else if (Super::i < Super::max_iter) {
> ret = false;
> diff --git a/itl/itl/preconditioner/ilu.hpp b/itl/itl/preconditioner/ilu.hpp
> index 9889614..b4e71e1 100644
> --- a/itl/itl/preconditioner/ilu.hpp
> +++ b/itl/itl/preconditioner/ilu.hpp
> @@ -6,8 +6,8 @@
> #ifndef ITL_ILU_HPP
> #define ITL_ILU_HPP
>
> -#include "itl/preconditioner/detail/triangular.hpp"
> #include "itl/preconditioner/detail/tri_solve.hpp"
> +#include "itl/preconditioner/detail/triangular.hpp"
>
> namespace itl {
>
--------------090106020405040407040902
Content-Type: text/plain; charset=UTF-8;
name="itl.h"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
filename="itl.h"
// -*- c++ -*-
//
//-----------------------------------------------------------------------
//
// $COPYRIGHT$
//
//=======================================================================
#ifndef ITL_ITL_H
#define ITL_ITL_H
/**@name Iterative Template Library
The following is the requirements for the parameter types
used in the ITL template functions.
The Matrix object must either follow the MTL requirements for a
matrix, be derived from the multiplier class (for matrix-free
multiplication), or have a specialized matvec::mult function defined
for it.
\begin{verbatim}
class Vector {
forward_iterator begin();
forward_iterator end();
};
\end{verbatim}
forward iterator refers to the iterator requirement
defined in the STL.
\begin{verbatim}
class Preconditioner {
void solve(const VectorX& x, VectorZ& z);
void trans_solve(const VectorX& x, VectorZ& z);
};
\end{verbatim}
The Preconditioner object performs a preconditioning operation based
on vector x and stores the result in vector z. The trans solve()
method only need be defined when the preconditioner is used with an
iterative solver that requires it.
\begin{verbatim}
class Iteration {
bool first();
bool finished(const VectorX& r);
bool finished(const Real& r);
bool converged(const VectorX& r);
bool converged(const Real& r);
void operator++();
void failed();
int error_code();
int iterations();
Real resid();
};
\end{verbatim}
The Iteration object calculates whether the solution has reached the
desired accuracy, or whether the maximum number of iterations has
been reached. The method finished() checks both convergence and
number of iterations. The method converged() only checks
convergence. The error code() method is used to determine the return
value for the this iterative solver function. The first() method is
used to determine the first iteration of the loop.
For all algorithms, if the error\_code() is 0, it suggests the algorithm
converges. Otherwise, if the error\_code() returns 1, it means the maximum
number of iteration has been reached but the desired accuacy is not reached.
For other return codes, see the respective document.
*/
/**
@name Utilities
@memo Utility classes and routines for ITL
*/
#include <itl/itl_config.h>
#include <iostream>
#include <complex>
#include <string>
namespace itl {
template <class Real>
class basic_iteration {
public:
typedef Real real;
template <class Vector>
basic_iteration(const Vector &b, int max_iter_, Real t, Real a = Real(0))
: error(0), i(0), normb_(std::abs(itl::two_norm(b))),
max_iter(max_iter_), rtol_(t), atol_(a) { }
basic_iteration(Real nb, int max_iter_, Real t, Real a = Real(0))
: error(0), i(0), normb_(nb), max_iter(max_iter_), rtol_(t), atol_(a) {}
template <class Vector>
bool finished(const Vector &r) {
Real normr_ = std::abs(itl::two_norm(r));
if (converged(normr_)) {
return true;
} else if (i < max_iter) {
return false;
} else {
error = 1;
return true;
}
}
bool finished(const Real &r) {
if (converged(r)) {
return true;
} else if (i < max_iter) {
return false;
} else {
error = 1;
return true;
}
}
template <typename T>
bool finished(const std::complex<T>& r) {
if (converged(std::abs(r))) {
return true;
} else if (i < max_iter) {
return false;
} else {
error = 1;
return true;
}
}
inline bool converged(const Real &r) {
resid_ = r / normb_;
return (resid_ <= rtol_ || r < atol_); //relative or absolute tolerance.
}
inline void operator++() { ++i; }
inline bool first() { return i == 0; }
inline int error_code() { return error; }
inline int iterations() { return i; }
inline Real resid() { return resid_ * normb_; }
inline Real normb() const { return normb_; }
inline Real tol() { return rtol_; }
inline Real atol() { return atol_; }
inline void fail(int err_code) { error = err_code; }
inline void fail(int err_code, const std::string &msg)
{ error = err_code; err_msg = msg; }
inline void set(Real v) { normb_ = v; }
protected:
int error;
int i;
const Real normb_;
int max_iter;
Real rtol_;
Real atol_;
Real resid_;
std::string err_msg;
};
template <class Real>
class noisy_iteration : public basic_iteration<Real> {
typedef basic_iteration<Real> Super;
int m_print_frequency;
bool m_verb;
void print_message(bool converged = false) {
if ( m_verb && ( !( Super::i % m_print_frequency ) || converged) )
std::cout << "Iteration " << Super::i << ": residual is "
<< Super::resid()
<< " \t relative error is " << Super::resid() / Super::normb()
<< std::endl;
}
public:
template <class Vector>
noisy_iteration(const Vector &b, int max_iter_,
Real tol_, Real atol_ = Real(0),
int print_frequency = 100,
bool verb = true)
: Super(b, max_iter_, tol_, atol_),
m_print_frequency(print_frequency),
m_verb(verb) {
if ( m_print_frequency == 0 ) {
m_print_frequency = 100;
}
}
noisy_iteration(Real nb, int max_iter_,
Real tol_, Real atol_ = Real(0),
int print_frequency = 100,
bool verb = true)
: Super(nb, max_iter_, tol_, atol_),
m_print_frequency(print_frequency),
m_verb(verb) {
if ( m_print_frequency == 0 ) {
m_print_frequency = 100;
}
}
template <class Vector>
bool finished(const Vector &r) {
using std::cout;
using std::endl;
Real normr_ = std::abs(itl::two_norm(r));
bool ret;
if (this->converged(normr_)) {
ret = true;
} else if (Super::i < Super::max_iter) {
ret = false;
} else {
Super::error = 1;
ret = true;
}
print_message(ret);
return ret;
}
bool finished(const Real &r) {
using std::cout;
using std::endl;
bool ret;
if (this->converged(r)) {
ret = true;
} else if (Super::i < Super::max_iter) {
ret = false;
} else {
Super::error = 1;
ret = true;
}
print_message(ret);
return ret;
}
template <typename T>
bool finished(const std::complex<T>& r) { //for the case of complex
using std::cout;
using std::endl;
bool ret;
if (this->converged(std::abs(r))) {
ret = true;
} else if (Super::i < Super::max_iter) {
ret = false;
} else {
Super::error = 1;
ret = true;
}
print_message(ret);
return ret;
}
int error_code() {
using std::cout;
using std::endl;
if (m_verb) {
if ( Super::error )
cout << "Warnning: linear solver stopped without required accuracy!"
<< " The error code = " << Super::error << endl;
else {
cout << "Linear solver has achieved the required accuracy. " << endl;
}
cout << Super::iterations() << " iterations" << endl;
cout << Super::resid() << " is actual final residual. " << endl
<< Super::resid() / Super::normb()
<< " is actual relative tolerance achieved. "
<< endl;
cout << "Relative tol: " << Super::rtol_
<< " Absolute tol: " << Super::atol_ << endl;
}
return Super::error;
}
};
struct identity_preconditioner {
identity_preconditioner operator()() const {
identity_preconditioner p;
return p;
}
identity_preconditioner left() const {
identity_preconditioner p;
return p;
}
identity_preconditioner right() const {
identity_preconditioner p;
return p;
}
};
template <class VecX, class VecZ>
inline void solve(const identity_preconditioner &M, const VecX &x,
const VecZ &z) {
itl::copy(x, const_cast<VecZ &>(z));
}
template <class VecX, class VecZ>
inline void trans_solve(const identity_preconditioner &M,
const VecX &x, const VecZ &z) {
itl::copy(x, const_cast<VecZ &>(z));
}
}
#endif
--------------090106020405040407040902
Content-Type: text/x-c++hdr;
name="ilu.hpp"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
filename="ilu.hpp"
//-----------------------------------------------------------------------
//
// $COPYRIGHT$
//
//=======================================================================
#ifndef ITL_ILU_HPP
#define ITL_ILU_HPP
#include "itl/preconditioner/detail/tri_solve.hpp"
#include "itl/preconditioner/detail/triangular.hpp"
namespace itl {
template <class T, class SizeT = size_t>
class ILU {
public:
typedef T value_type;
typedef SizeT size_type;
private:
std::vector<value_type> L_val;
std::vector<size_type> L_ind;
std::vector<size_type> L_ptr;
std::vector<value_type> U_val;
std::vector<size_type> U_ind;
std::vector<size_type> U_ptr;
size_type n;
public:
template <class Matrix>
ILU(const Matrix &A) {
itl::matrix_split(A, L_val, L_ind, L_ptr,
U_val, U_ind, U_ptr);
n = L_ptr.size() - 1;
size_type i, j, qn, pn, rn;
for (i = 1; i < n; i++) {
for (j = L_ptr[i]; j < L_ptr[i + 1]; j++) {
pn = U_ptr[L_ind[j]];
value_type multiplier = (L_val[j] /= U_val[pn]);
qn = j + 1;
rn = U_ptr[i];
for (pn++; U_ind[pn] < i && pn < U_ptr[L_ind[j] + 1]; pn++) {
while (L_ind[qn] < U_ind[pn] && qn < L_ptr[i + 1]) {
qn++;
}
if (U_ind[pn] == L_ind[qn] && qn < L_ptr[i + 1]) {
L_val[qn] -= multiplier * U_val[pn];
}
}
for (; pn < U_ptr[L_ind[j] + 1]; pn++) {
while (U_ind[rn] < U_ind[pn] && rn < U_ptr[i + 1]) {
rn++;
}
if (U_ind[pn] == U_ind[rn] && rn < U_ptr[i + 1]) {
U_val[rn] -= multiplier * U_val[pn];
}
}
}
}
}
void print() {
for (size_type i = 0; i < n; ++i) {
size_type b = L_ptr[i], e = L_ptr[i + 1];
for (size_type j = b; j < e; ++j) {
std::cout << L_val[j] << "[" << i << "," << L_ind[j] << "] ";
}
std::cout << "|| ";
b = U_ptr[i], e = U_ptr[i + 1];
for (size_type j = b; j < e; ++j) {
std::cout << U_val[j] << "[" << i << "," << U_ind[j] << "] ";
}
std::cout << std::endl;
}
}
template <class U>
void solve(U *x) const {
{
bool Unit = true;
bool Upper = false;
bool Trans = false;
triangular_solve(n, &L_val[0], &L_ind[0], &L_ptr[0], x,
Unit, Upper, Trans);
}
{
bool Unit = false;
bool Upper = true;
bool Trans = false;
triangular_solve(n, &U_val[0], &U_ind[0], &U_ptr[0], x,
Unit, Upper, Trans);
}
}
template <class U>
void trans_solve(U *x) const {
{
bool Unit = false;
bool Upper = true;
bool Trans = true;
triangular_solve(n, &U_val[0], &U_ind[0], &U_ptr[0], x,
Unit, Upper, Trans);
}
{
bool Unit = true;
bool Upper = false;
bool Trans = true;
triangular_solve(n, &L_val[0], &L_ind[0], &L_ptr[0], x,
Unit, Upper, Trans);
}
}
Triangular<value_type, size_type> left() const
{ return Triangular<value_type, size_type>(L_val, L_ind, L_ptr, true, false, false); }
Triangular<value_type, size_type> right() const
{ return Triangular<value_type, size_type>(U_val, U_ind, U_ptr, false, true, false); }
};
template <class T, class S, class VectorB, class VectorX>
inline void solve(const ILU<T, S>& M, const VectorB &b, const VectorX &x)
{ itl::copy(b, const_cast<VectorX &>(x)); M.solve(get_data(x)); }
template <class T, class S, class VectorB, class VectorX>
inline void trans_solve(const ILU<T, S>& M, const VectorB &b,
const VectorX &x)
{ itl::copy(b, const_cast<VectorX &>(x)); M.trans_solve(get_data(x)); }
}
#endif
--------------090106020405040407040902
Content-Type: text/plain; charset="us-ascii"
MIME-Version: 1.0
Content-Transfer-Encoding: 7bit
Content-Disposition: inline
_______________________________________________
This list is archived at http://www.osl.iu.edu/MailArchives/mtl-devel/
--------------090106020405040407040902--