Skip to content
Snippets Groups Projects
Commit a88d34f3 authored by Uli Sack's avatar Uli Sack Committed by usack
Browse files

introduce implementation of normSquared to avoid often unnecessary roots and...

introduce implementation of normSquared to avoid often unnecessary roots and squaring;  compute norm by sqrt(normSquared); double->field_type; harmonized argument names with base class (i.e. v->f)

[[Imported from SVN: r12076]]
parent 59281200
Branches
Tags
No related merge requests found
...@@ -22,31 +22,39 @@ template <class LowRankFactor=Dune::BlockVector<Dune::FieldVector <double,1> >, ...@@ -22,31 +22,39 @@ template <class LowRankFactor=Dune::BlockVector<Dune::FieldVector <double,1> >,
class FullNorm: public Norm<VectorType> class FullNorm: public Norm<VectorType>
{ {
public: public:
FullNorm(const double alpha, const LowRankFactor &lowRankFactor) : typedef typename VectorType::field_type field_type;
FullNorm(const field_type alpha, const LowRankFactor &lowRankFactor) :
lowRankFactor_(lowRankFactor), lowRankFactor_(lowRankFactor),
alpha(alpha) alpha(alpha)
{} {}
//! Compute the norm of the given vector //! Compute the norm of the given vector
double operator()(const VectorType &v) const field_type operator()(const VectorType &f) const
{
return std::sqrt(this->FullNorm<LowRankFactor,VectorType>::normSquared(f));
}
//! Compute the square of the norm of the given vector
virtual field_type normSquared(const VectorType& f) const
{ {
VectorType r(lowRankFactor_.N()); VectorType r(lowRankFactor_.N());
r = 0.0; r = 0.0;
lowRankFactor_.umv(v,r); lowRankFactor_.umv(f,r);
return std::sqrt(std::abs(alpha*(r*r))); return std::abs(alpha*(r*r));
} }
//! Compute the norm of the difference of two vectors //! Compute the norm of the difference of two vectors
double diff(const VectorType &v1, const VectorType &v2) const field_type diff(const VectorType &f1, const VectorType &f2) const
{ {
VectorType r(lowRankFactor_.N()); VectorType r(lowRankFactor_.N());
r = 0.0; r = 0.0;
for (typename LowRankFactor::size_type k=0; k<lowRankFactor_.N(); ++k) for (typename LowRankFactor::size_type k=0; k<lowRankFactor_.N(); ++k)
for (typename LowRankFactor::size_type i=0; i<lowRankFactor_.M(); ++i) for (typename LowRankFactor::size_type i=0; i<lowRankFactor_.M(); ++i)
lowRankFactor_[k][i].umv(v1[i]-v2[i],r[k]); lowRankFactor_[k][i].umv(f1[i]-f2[i],r[k]);
return std::sqrt(std::abs(alpha*(r*r))); return std::sqrt(std::abs(alpha*(r*r)));
} }
...@@ -54,7 +62,7 @@ class FullNorm: public Norm<VectorType> ...@@ -54,7 +62,7 @@ class FullNorm: public Norm<VectorType>
private: private:
const LowRankFactor& lowRankFactor_; const LowRankFactor& lowRankFactor_;
const double alpha; const field_type alpha;
}; };
...@@ -65,30 +73,38 @@ class FullNorm<Dune::BlockVector<Dune::FieldVector<double,1> >, Dune::BlockVecto ...@@ -65,30 +73,38 @@ class FullNorm<Dune::BlockVector<Dune::FieldVector<double,1> >, Dune::BlockVecto
typedef Dune::BlockVector<Dune::FieldVector<double,1> > VectorType; typedef Dune::BlockVector<Dune::FieldVector<double,1> > VectorType;
typedef VectorType::size_type SizeType; typedef VectorType::size_type SizeType;
typedef double field_type;
public: public:
FullNorm(const double alpha, const VectorType &m) : FullNorm(const field_type alpha, const VectorType &m) :
m(m), m(m),
alpha(alpha) alpha(alpha)
{} {}
//! Compute the norm of the given vector //! Compute the norm of the given vector
double operator()(const VectorType &v) const field_type operator()(const VectorType &f) const
{
return std::sqrt(this->FullNorm<>::normSquared(f));
}
//! Compute the square of the norm of the given vector
virtual field_type normSquared(const VectorType& f) const
{ {
double r = 0.0; field_type r = 0.0;
for (SizeType row = 0; row < v.size(); ++row) for (SizeType row = 0; row < f.size(); ++row)
r += m[row] * v[row]; r += m[row] * f[row];
return std::sqrt(std::abs(alpha*r*r)); return std::abs(alpha*r*r);
} }
//! Compute the norm of the difference of two vectors //! Compute the norm of the difference of two vectors
double diff(const VectorType &v1, const VectorType &v2) const field_type diff(const VectorType &f1, const VectorType &f2) const
{ {
double r = 0.0; field_type r = 0.0;
for (SizeType row = 0; row < v1.size(); ++row) for (SizeType row = 0; row < f1.size(); ++row)
r += (double)m[row] * (v1[row] - v2[row]); r += (field_type)m[row] * (f1[row] - f2[row]);
return std::sqrt(std::abs(alpha*r*r)); return std::sqrt(std::abs(alpha*r*r));
} }
...@@ -96,7 +112,7 @@ class FullNorm<Dune::BlockVector<Dune::FieldVector<double,1> >, Dune::BlockVecto ...@@ -96,7 +112,7 @@ class FullNorm<Dune::BlockVector<Dune::FieldVector<double,1> >, Dune::BlockVecto
private: private:
const VectorType &m; const VectorType &m;
const double alpha; const field_type alpha;
}; };
#endif #endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment