source: trunk/yat/utility/vector.cc @ 1027

Last change on this file since 1027 was 1027, checked in by Peter, 14 years ago

going back to previous design in which view and const_view are in different classes. Having them in same didnt work as expected. There is a problem converting vector::iterator to vector::const_iterator. I'll open a separate ticket for this issue.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 6.0 KB
Line 
1// $Id: vector.cc 1027 2008-02-02 21:29:29Z peter $
2
3/*
4  Copyright (C) 2003 Daniel Dalevi, Peter Johansson
5  Copyright (C) 2004 Jari Häkkinen, Peter Johansson
6  Copyright (C) 2005, 2006 Jari Häkkinen, Markus Ringnér, Peter Johansson
7  Copyright (C) 2007 Jari Häkkinen, Peter Johansson
8
9  This file is part of the yat library, http://trac.thep.lu.se/yat
10
11  The yat library is free software; you can redistribute it and/or
12  modify it under the terms of the GNU General Public License as
13  published by the Free Software Foundation; either version 2 of the
14  License, or (at your option) any later version.
15
16  The yat library is distributed in the hope that it will be useful,
17  but WITHOUT ANY WARRANTY; without even the implied warranty of
18  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19  General Public License for more details.
20
21  You should have received a copy of the GNU General Public License
22  along with this program; if not, write to the Free Software
23  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
24  02111-1307, USA.
25*/
26
27#include "vector.h"
28#include "matrix.h"
29#include "utility.h"
30#include "yat/random/random.h"
31
32#include <algorithm>
33#include <cassert>
34#include <cmath>
35#include <iostream>
36#include <sstream>
37#include <utility>
38#include <vector>
39
40namespace theplu {
41namespace yat {
42namespace utility {
43
44
45  vector::vector(void)
46    : VectorMutable()
47  {
48  }
49
50
51  vector::vector(size_t n, double init_value)
52    : VectorMutable(gsl_vector_alloc(n))
53  {
54    if (!vec_)
55      throw utility::GSL_error("vector::vector failed to allocate memory");
56    assert(const_vec_);
57    all(init_value);
58  }
59
60
61  vector::vector(const vector& other)
62    : VectorMutable(create_gsl_vector_copy(other))
63  {
64  }
65
66
67  vector::vector(const VectorBase& other)
68    : VectorMutable(create_gsl_vector_copy(other))
69  {
70  }
71
72
73  vector::vector(std::istream& is, char sep)
74    throw (utility::IO_error, std::exception)
75    : VectorMutable()
76  {
77    // read the data file and store in stl vectors (dynamically
78    // expandable)
79    std::vector<std::vector<double> > data_matrix;
80    u_int nof_columns=0;
81    u_int nof_rows=0;
82    std::string line;
83    while(getline(is, line, '\n')) {
84      // Empty lines
85      if (!line.size())
86          continue;
87      nof_rows++;
88
89      std::vector<double> v;
90      std::string element;
91      std::stringstream ss(line);
92      bool ok=true;
93      while(ok) {
94        if(sep=='\0')
95          ok=(ss>>element);
96        else 
97          ok=getline(ss, element, sep);
98        if(!ok)
99          break;       
100
101        if(utility::is_double(element)) {
102          v.push_back(atof(element.c_str()));
103        }
104        else if (!element.size() || utility::is_nan(element)) {
105          v.push_back(std::numeric_limits<double>::quiet_NaN());
106        }
107        else {
108          std::stringstream ss("Warning: '");
109          ss << element << "' is not a double.";
110          throw IO_error(ss.str());
111        }
112      } 
113      if(sep!='\0' && line[line.size()-1]==sep) // add NaN for final separator
114          v.push_back(std::numeric_limits<double>::quiet_NaN());
115      if (!nof_columns)
116        nof_columns=v.size();
117      else if (nof_rows && (nof_columns>1)) {
118        std::ostringstream s;
119        s << "vector::vector(std::istream&) data file error:\n"
120          << "    File has inconsistent number of rows (" << nof_rows
121          << ") and columns (" << nof_columns
122          << ").\n    Expected a row or a column vector.";
123        throw utility::IO_error(s.str());
124      }
125      else if (v.size()!=nof_columns) {
126        std::ostringstream s;
127        s << "vector::vector(std::istream&) data file error:\n"
128          << "    Line " << nof_rows << " has " << v.size()
129          << " columns; expected " << nof_columns << " column.";
130        throw utility::IO_error(s.str());
131      }
132      data_matrix.push_back(v);
133    }
134
135    // manipulate the state of the stream to be good
136    is.clear(std::ios::goodbit);
137    // convert the data to a gsl vector
138    vec_ = gsl_vector_alloc(nof_rows*nof_columns);
139    if (!vec_)
140      throw utility::GSL_error("vector::vector failed to allocate memory");
141    size_t n=0;
142    // if gsl error handler disabled, out of bounds index will not
143    // abort the program.
144    for (size_t i=0; i<nof_rows; i++)
145      for (size_t j=0; j<nof_columns; j++) 
146        gsl_vector_set( vec_, n++, data_matrix[i][j] );
147    const_vec_ = vec_;
148  }
149
150
151  vector::~vector(void)
152  {
153    delete_allocated_memory();
154  }
155
156
157  gsl_vector* vector::create_gsl_vector_copy(const VectorBase& other) const
158  {
159    gsl_vector* vec = gsl_vector_alloc(other.size());
160    if (!vec)
161      throw utility::GSL_error("vector::create_gsl_vector_copy failed to allocate memory");
162    if (gsl_vector_memcpy(vec, other.gsl_vector_p()))
163      throw utility::GSL_error("vector::create_gsl_matrix_copy dimension mis-match");
164    return vec;
165  }
166
167
168  void vector::delete_allocated_memory(void) 
169  {
170    if (vec_)
171      gsl_vector_free(vec_);
172    const_vec_ = vec_ = NULL;
173  }
174
175
176  bool vector::isview(void) const
177  {
178    return false;
179  }
180
181
182  void vector::resize(size_t n, double init_value)
183  {
184    delete_allocated_memory();
185    const_vec_ = vec_ = gsl_vector_alloc(n);
186    if (!vec_)
187      throw utility::GSL_error("vector::vector failed to allocate memory");
188    all(init_value);
189   
190  } 
191
192
193  /*
194  //Peter use swap idiom
195  const VectorBase& vector::operator=( const VectorBase& other )
196  {
197    if (!other.size())
198      vec_=NULL;
199    else  {   
200      if (size()!=other.size())
201        resize(other.size(),0.0);
202      gsl_vector_memcpy(vec_,other.gsl_vector_p());
203    }
204    return *this;
205  }
206  */
207
208  const vector& vector::operator=( const vector& other )
209  {
210    return assign(other);
211  }
212
213  /* 
214  const VectorBase& vector::assign(VectorBase& other)
215  {
216    assign(other);
217  }
218  */
219
220  const vector& vector::assign(const VectorBase& other)
221  {
222    if (!other.size())
223      vec_=NULL;
224    else  {   
225      if (size()!=other.size())
226        resize(other.size(),0.0);
227      gsl_vector_memcpy(vec_,other.gsl_vector_p());
228    }
229    const_vec_ = vec_;
230    return *this;
231  } 
232
233
234  void set_basis(vector& v, size_t i)
235  {
236    assert(v.gsl_vector_p());
237    gsl_vector_set_basis(v.gsl_vector_p(),i);
238  }
239
240
241  void sort(vector& v)
242  {
243    assert(v.gsl_vector_p());
244    std::sort(v.begin(), v.end());
245  }
246
247
248  void swap(vector& v, vector& w)
249  {
250    assert(v.gsl_vector_p()); assert(w.gsl_vector_p());
251    int status=gsl_vector_swap(v.gsl_vector_p(),w.gsl_vector_p());
252    if (status)
253      throw utility::GSL_error(std::string("swap(vector&,vector&)",status));
254  }
255
256
257}}} // of namespace utility, yat, and thep
Note: See TracBrowser for help on using the repository browser.