source: trunk/src/vector.cc @ 227

Last change on this file since 227 was 227, checked in by Jari Häkkinen, 17 years ago

Started reimplementation of the vector class.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.1 KB
Line 
1// $Id: vector.cc 227 2005-02-01 00:52:27Z jari $
2
3#include "vector.h"
4#include "random_singleton.h"
5
6#include <iostream>
7#include <sstream>
8#include <string>
9#include <vector>
10#include <utility>
11
12
13
14
15
16namespace theplu {
17namespace gslapi {
18
19
20
21  vector::vector(void)
22    : v_(NULL), view_(NULL)
23  {
24  }
25
26
27
28  vector::vector(const size_t n,const double init_value)
29    : view_(NULL)
30  {
31    v_ = gsl_vector_alloc(n);
32    set_all(init_value);
33  }
34
35
36
37  vector::vector(const vector& other)
38    : view_(NULL)
39  {
40    v_ = other.TEMP_gsl_vector_copy();
41  }
42
43
44
45  vector::vector(gsl_vector* vec)
46    : v_(vec), view_(NULL)
47  {
48  }
49
50
51
52  vector::vector(std::istream& is)
53    : view_(NULL)
54  {
55    using namespace std;
56 
57    // read the data file and store in stl vectors (dynamically expandable)
58    std::vector<std::vector<double> > data_matrix;
59    u_int nof_columns=0;
60    u_int nof_rows = 0;
61    string s;
62    for (nof_rows = 0; getline(is, s, '\n'); nof_rows++) {
63      istringstream line(s);
64      std::vector<double> v;
65      string tmp_string;
66      while (line >> tmp_string) {
67        if(!is.good()) {
68          // Jari. change to throw exception, remove unnecessary includes
69          cerr << "vector::vector(std::istream&): "
70               << "error reading data file!" << endl;
71          exit(1);
72        }
73        double t=atof(tmp_string.c_str());
74        // Here we should check that it actually was a double!!!!!
75        v.push_back(t);
76      }
77
78      // Ignoring empty lines
79      if(!v.size()){
80        nof_rows--;
81        continue;
82      }
83
84      if(nof_columns==0)
85        nof_columns=v.size();
86      else if(v.size()!=nof_columns) {
87        // Jari. change to throw exception, remove unnecessary includes
88        cerr << "vector::vector(std::istream&) data file error: "
89             << "line" << nof_rows+1 << " has " << v.size() 
90             << " columns; expected " << nof_columns
91             << " columns"
92             << endl; 
93        exit(1);
94      }
95      data_matrix.push_back(v);
96    }
97 
98    // manipulate the state of the stream to be good
99    is.clear(std::ios::goodbit);
100 
101    // convert the data to a gsl vector and check that data file is a
102    // column vector or a row vector
103    if(nof_columns==1) {
104      v_ = gsl_vector_alloc ( nof_rows );
105      for(u_int i=0;i<nof_rows;i++) 
106        gsl_vector_set( v_, i, data_matrix[i][0] );
107    }
108    else if(nof_rows==1){
109      v_ = gsl_vector_alloc ( nof_columns );
110      for(u_int i=0;i<nof_columns;i++) 
111        gsl_vector_set( v_, i, data_matrix[0][i] );
112    }
113    else {
114      // Jari. change to throw exception, remove unnecessary includes
115      cerr << "vector::vector(std::istream&) data file error: "
116           << "file has " << nof_rows << " rows and " << nof_columns
117           << " columns; expected a row vector or a column vector"
118           << endl;
119      exit(1);
120    }
121  }
122
123
124
125  vector::~vector(void)
126  {
127    if (v_) {
128      gsl_vector_free(v_);
129      v_=NULL;
130    }
131    delete view_;
132  }
133
134
135
136  gsl_vector* vector::TEMP_gsl_vector_copy(void) const
137  {
138    gsl_vector* vec = gsl_vector_alloc(size());
139    gsl_vector_memcpy(vec,v_);
140    return vec;
141  }
142
143
144
145  std::pair<double,double> vector::minmax(void) const
146  {
147    double min, max;
148    gsl_vector_minmax(v_,&min,&max);
149    return std::pair<double,double>(min,max);
150  }
151
152
153
154  std::pair<size_t,size_t> vector::minmax_index(void) const
155  {
156    size_t min_index, max_index;
157    gsl_vector_minmax_index(v_,&min_index,&max_index);
158    return std::pair<size_t,size_t>(min_index,max_index);
159  }
160
161
162
163/* Jari, remove this section
164  std::pair<size_t,size_t> vector::minmax_index(const std::vector<size_t>& subset ) const
165  {
166    size_t min_index = subset[0];
167    size_t max_index = subset[0];
168    for (unsigned int i=1; i<subset.size(); i++) {
169      if (gsl_vector_get( v_, subset[i]) < gsl_vector_get( v_, min_index))
170        min_index = subset[i];
171      else if (gsl_vector_get( v_, subset[i]) > gsl_vector_get( v_, max_index))
172        max_index = subset[i];
173    }
174    return std::pair<size_t,size_t>(min_index, max_index);
175  }
176*/
177
178  void vector::shuffle(void) const
179  {
180    vector vec(*this);
181    theplu::cpptools::random_singleton* 
182      rnd=theplu::cpptools::random_singleton::get_instance();
183    for (size_t i=0; i<vec.size(); i++){
184      size_t index = rnd->get_uniform_int(vec.size()-i);
185      gsl_vector_set(v_,i,vec(index));
186      vec(index)=vec(vec.size()-i-1);
187    }
188  }
189
190  double vector::sum(void) const
191  {
192    double sum = 0;
193    for (size_t i=0; i<size(); i++)
194      sum += gsl_vector_get( v_, i );
195    return( sum );
196  } 
197
198
199
200  vector vector::operator-(void) const
201  {
202    vector res( *this );
203    gsl_vector_scale (res.v_,-1.);
204    return res;
205  }
206
207
208
209  double vector::operator*( const vector &other ) const
210  {
211    // Jari, check for gsl_support
212    double res = 0.0;;
213    for ( size_t i = 0; i < size(); ++i ) 
214      res += other[i] * (*this)[i];
215    return res;
216  }
217
218
219
220  vector vector::operator+(const vector &other) const
221  {
222    vector res(*this);
223    gsl_vector_add(res.v_,other.v_);
224    return res;
225  }
226
227
228
229  vector vector::operator-( const vector &other ) const
230  {
231    vector res( *this );
232    gsl_vector_sub(res.v_,other.v_);
233    return res;
234  }
235
236
237
238  vector& vector::operator=( const vector& other )
239  {
240    if( this != &other ) {
241      if ( v_ )
242        gsl_vector_free( v_ );
243      v_ = other.TEMP_gsl_vector_copy();
244    }
245    return *this;
246  } 
247
248
249
250  std::ostream& theplu::gslapi::operator<<(std::ostream& s, const vector& a)
251  {
252    s.setf(std::ios::fixed);
253
254    for (size_t j = 0; j < a.size(); ++j) {
255      s << a[j];
256      if ( (j+1)<a.size() )
257        s << " ";
258    }
259
260    return s;
261  }
262
263
264}} // of namespace gslapi and namespace thep
Note: See TracBrowser for help on using the repository browser.