← Back to branch summary

~ubuntu-branches/ubuntu/intrepid/libmesh/intrepid

« back to all changes in this revision

Viewing changes to src/numerics/numeric_vector.C

  • Committer: Bazaar Package Importer
  • Author(s): Christophe Prud'homme
  • Date: 2007-09-03 22:11:32 UTC
  • Revision ID: james.westby@ubuntu.com-20070903221132-9debourfspc47qwt
Tags: upstream-0.6.0~rc2.dfsg
Import upstream version 0.6.0~rc2.dfsg

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/numerics/numeric_vector.C
 
1
// $Id: numeric_vector.C,v 1.18 2005/11/30 00:28:32 roystgnr Exp $
 
2
 
 
3
// The libMesh Finite Element Library.
 
4
// Copyright (C) 2002-2005  Benjamin S. Kirk, John W. Peterson
 
5
  
 
6
// This library is free software; you can redistribute it and/or
 
7
// modify it under the terms of the GNU Lesser General Public
 
8
// License as published by the Free Software Foundation; either
 
9
// version 2.1 of the License, or (at your option) any later version.
 
10
  
 
11
// This library is distributed in the hope that it will be useful,
 
12
// but WITHOUT ANY WARRANTY; without even the implied warranty of
 
13
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 
14
// Lesser General Public License for more details.
 
15
  
 
16
// You should have received a copy of the GNU Lesser General Public
 
17
// License along with this library; if not, write to the Free Software
 
18
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 
19
 
 
20
 
 
21
 
 
22
// C++ includes
 
23
#include <cmath> // for std::abs
 
24
 
 
25
// Local Includes
 
26
#include "numeric_vector.h"
 
27
#include "laspack_vector.h"
 
28
#include "petsc_vector.h"
 
29
 
 
30
 
 
31
 
 
32
//------------------------------------------------------------------
 
33
// NumericVector methods
 
34
 
 
35
// Full specialization for Real datatypes
 
36
template <typename T>
 
37
AutoPtr<NumericVector<T> >
 
38
NumericVector<T>::build(const SolverPackage solver_package)
 
39
{
 
40
  // Build the appropriate vector
 
41
  switch (solver_package)
 
42
    {
 
43
 
 
44
 
 
45
#ifdef HAVE_LASPACK
 
46
    case LASPACK_SOLVERS:
 
47
      {
 
48
        AutoPtr<NumericVector<T> > ap(new LaspackVector<T>);
 
49
        return ap;
 
50
      }
 
51
#endif
 
52
 
 
53
 
 
54
#ifdef HAVE_PETSC
 
55
    case PETSC_SOLVERS:
 
56
      {
 
57
        AutoPtr<NumericVector<T> > ap(new PetscVector<T>);
 
58
        return ap;
 
59
      }
 
60
#endif
 
61
 
 
62
    default:
 
63
      std::cerr << "ERROR:  Unrecognized solver package: "
 
64
                << solver_package
 
65
                << std::endl;
 
66
      error();
 
67
    }
 
68
    
 
69
  AutoPtr<NumericVector<T> > ap(NULL);
 
70
  return ap;    
 
71
}
 
72
 
 
73
 
 
74
 
 
75
// Full specialization for float datatypes (DistributedVector wants this)
 
76
 
 
77
template <>
 
78
int NumericVector<float>::compare (const NumericVector<float> &other_vector,
 
79
                                   const Real threshold) const
 
80
{
 
81
  assert (this->initialized());
 
82
  assert (other_vector.initialized());
 
83
  assert (this->first_local_index() == other_vector.first_local_index());
 
84
  assert (this->last_local_index()  == other_vector.last_local_index());
 
85
 
 
86
  int rvalue     = -1;
 
87
  unsigned int i = first_local_index();
 
88
 
 
89
  do
 
90
    {
 
91
      if ( std::abs( (*this)(i) - other_vector(i) ) > threshold )
 
92
        rvalue = i;
 
93
      else
 
94
        i++;
 
95
    }
 
96
  while (rvalue==-1 && i<last_local_index());
 
97
 
 
98
  return rvalue;
 
99
}
 
100
 
 
101
// Full specialization for double datatypes
 
102
template <>
 
103
int NumericVector<double>::compare (const NumericVector<double> &other_vector,
 
104
                                    const Real threshold) const
 
105
{
 
106
  assert (this->initialized());
 
107
  assert (other_vector.initialized());
 
108
  assert (this->first_local_index() == other_vector.first_local_index());
 
109
  assert (this->last_local_index()  == other_vector.last_local_index());
 
110
 
 
111
  int rvalue     = -1;
 
112
  unsigned int i = first_local_index();
 
113
 
 
114
  do
 
115
    {
 
116
      if ( std::abs( (*this)(i) - other_vector(i) ) > threshold )
 
117
        rvalue = i;
 
118
      else
 
119
        i++;
 
120
    }
 
121
  while (rvalue==-1 && i<last_local_index());
 
122
 
 
123
  return rvalue;
 
124
}
 
125
 
 
126
#ifdef TRIPLE_PRECISION
 
127
// Full specialization for long double datatypes
 
128
template <>
 
129
int NumericVector<long double>::compare (const NumericVector<long double> &other_vector,
 
130
                                         const Real threshold) const
 
131
{
 
132
  assert (this->initialized());
 
133
  assert (other_vector.initialized());
 
134
  assert (this->first_local_index() == other_vector.first_local_index());
 
135
  assert (this->last_local_index()  == other_vector.last_local_index());
 
136
 
 
137
  int rvalue     = -1;
 
138
  unsigned int i = first_local_index();
 
139
 
 
140
  do
 
141
    {
 
142
      if ( std::abs( (*this)(i) - other_vector(i) ) > threshold )
 
143
        rvalue = i;
 
144
      else
 
145
        i++;
 
146
    }
 
147
  while (rvalue==-1 && i<last_local_index());
 
148
 
 
149
  return rvalue;
 
150
}
 
151
#endif
 
152
 
 
153
 
 
154
// Full specialization for Complex datatypes
 
155
template <>
 
156
int NumericVector<Complex>::compare (const NumericVector<Complex> &other_vector,
 
157
                                     const Real threshold) const
 
158
{
 
159
  assert (this->initialized());
 
160
  assert (other_vector.initialized());
 
161
  assert (this->first_local_index() == other_vector.first_local_index());
 
162
  assert (this->last_local_index()  == other_vector.last_local_index());
 
163
 
 
164
  int rvalue     = -1;
 
165
  unsigned int i = first_local_index();
 
166
 
 
167
  do
 
168
    {
 
169
      if (( std::abs( (*this)(i).real() - other_vector(i).real() ) > threshold ) ||
 
170
          ( std::abs( (*this)(i).imag() - other_vector(i).imag() ) > threshold ))
 
171
        rvalue = i;
 
172
      else
 
173
        i++;
 
174
    }
 
175
  while (rvalue==-1 && i<this->last_local_index());
 
176
 
 
177
  return rvalue;
 
178
}
 
179
 
 
180
 
 
181
 
 
182
//------------------------------------------------------------------
 
183
// Explicit instantiations
 
184
template class NumericVector<Number>;