source: trunk/yat/utility/Vector.cc @ 1929

Last change on this file since 1929 was 1929, checked in by Peter, 13 years ago

refs #510. using load function in Vector(istream) constructor

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.5 KB
Line 
1// $Id: Vector.cc 1929 2009-04-30 17:39:06Z 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, 2007 Jari Häkkinen, Peter Johansson, Markus Ringnér
7  Copyright (C) 2008 Jari Häkkinen, Peter Johansson
8  Copyright (C) 2009 Peter Johansson
9
10  This file is part of the yat library, http://dev.thep.lu.se/yat
11
12  The yat library is free software; you can redistribute it and/or
13  modify it under the terms of the GNU General Public License as
14  published by the Free Software Foundation; either version 3 of the
15  License, or (at your option) any later version.
16
17  The yat library is distributed in the hope that it will be useful,
18  but WITHOUT ANY WARRANTY; without even the implied warranty of
19  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20  General Public License for more details.
21
22  You should have received a copy of the GNU General Public License
23  along with yat. If not, see <http://www.gnu.org/licenses/>.
24*/
25
26#include "Vector.h"
27#include "utility.h"
28
29#include <algorithm>
30#include <cassert>
31#include <cmath>
32#include <iostream>
33#include <string>
34#include <sstream>
35#include <utility>
36#include <vector>
37
38namespace theplu {
39namespace yat {
40namespace utility {
41
42
43  Vector::Vector(void)
44    : VectorMutable()
45  {
46  }
47
48
49  Vector::Vector(size_t n, double init_value)
50    : VectorMutable(create_gsl_vector(n))
51  {
52    if (n)
53      all(init_value);
54    assert(n==0 || vec_);
55    assert(n==0 || const_vec_);
56  }
57
58
59  Vector::Vector(const Vector& other)
60    : VectorMutable(create_gsl_vector_copy(other))
61  {
62  }
63
64
65  Vector::Vector(const VectorBase& other)
66    : VectorMutable(create_gsl_vector_copy(other))
67  {
68  }
69
70
71  Vector::Vector(std::istream& is, char sep)
72    throw (utility::IO_error, std::exception)
73    : VectorMutable()
74  {
75    if (!is.good())
76      throw utility::IO_error("Vector: istream is not good");
77
78    // read the data file and store in stl vectors (dynamically
79    // expandable)
80    std::vector<std::vector<double> > data_matrix;
81    try {
82      load(is, data_matrix, sep, '\n', true);
83    }
84    catch (utility::IO_error& e) {
85      std::stringstream ss(e.what());
86      ss << "\nVector(std::istream&): invalid dimensions\n";
87      throw IO_error(ss.str());
88    }
89    catch (std::runtime_error& e) {
90      std::stringstream ss(e.what());
91      ss << "\nVector(std::istream&): invalid vector element\n";
92      throw IO_error(ss.str());
93    }
94
95    unsigned int nof_rows = data_matrix.size();
96    // if stream was empty, create nothing
97    if (!nof_rows)
98      return;
99
100    unsigned int nof_columns=data_matrix[0].size();
101
102    if (!(nof_rows==1 || nof_columns==1)) {
103      std::stringstream ss;
104      ss << "\nVector(std::istream&) data file error:\n"
105         << "    File has inconsistent number of rows (" << nof_rows
106         << ") and columns (" << nof_columns
107         << ").\n    Expected a row or column Vector.";
108      throw IO_error(ss.str());
109    }
110
111    // convert the data to a gsl vector
112    vec_ = gsl_vector_alloc(nof_rows*nof_columns);
113    if (!vec_)
114      throw utility::GSL_error("Vector::Vector failed to allocate memory");
115    size_t n=0;
116    // if gsl error handler disabled, out of bounds index will not
117    // abort the program.
118    for (size_t i=0; i<nof_rows; i++)
119      for (size_t j=0; j<nof_columns; j++) 
120        gsl_vector_set( vec_, n++, data_matrix[i][j] );
121    const_vec_ = vec_;
122  }
123
124
125  Vector::~Vector(void)
126  {
127    delete_allocated_memory();
128  }
129
130
131  const Vector& Vector::assign(const VectorBase& other)
132  {
133    // avoid self assignment
134    if (this == &other)
135      return *this;
136    if (!other.size()){
137      delete_allocated_memory();
138      return *this;
139    }
140    // indirect self assignment if begin <= other.begin() < end
141    if (begin() <= other.begin() && other.begin()<end()) {
142      if (size() != other.size()) {
143        gsl_vector* tmp = create_gsl_vector_copy(other);
144        delete_allocated_memory();
145        const_vec_ = vec_ = tmp;
146      }
147      return *this;
148    } 
149    if (size()==other.size()){
150      if (gsl_vector_memcpy(vec_, other.gsl_vector_p()))
151        throw utility::GSL_error("Vector::assign memcpy failed");
152    }
153    else {
154      delete_allocated_memory();
155      vec_ = create_gsl_vector_copy(other);
156    }
157    const_vec_ = vec_;
158    return *this;
159  } 
160
161
162  gsl_vector* Vector::create_gsl_vector_copy(const VectorBase& other) const
163  {
164    if (!other.size())
165      return NULL;
166    gsl_vector* vec = gsl_vector_alloc(other.size());
167    if (!vec)
168      throw utility::GSL_error("Vector::create_gsl_vector_copy failed to allocate memory");
169    if (gsl_vector_memcpy(vec, other.gsl_vector_p()))
170      throw utility::GSL_error("Vector::create_gsl_vector_copy memcpy failed");
171    return vec;
172  }
173
174
175  gsl_vector* Vector::create_gsl_vector(size_t n) const
176  {
177    if (!n)
178      return NULL;
179    gsl_vector* vec = gsl_vector_alloc(n);
180    if (!vec)
181      throw utility::GSL_error("Vector::create_gsl_vector failed to allocate memory");
182    return vec;
183  }
184
185
186  void Vector::delete_allocated_memory(void) 
187  {
188    if (vec_)
189      gsl_vector_free(vec_);
190    const_vec_ = vec_ = NULL;
191  }
192
193
194  bool Vector::isview(void) const
195  {
196    return false;
197  }
198
199
200  void Vector::resize(size_t n, double init_value)
201  {
202    delete_allocated_memory();
203    if (!n)
204      return;
205    const_vec_ = vec_ = gsl_vector_alloc(n);
206    if (!vec_)
207      throw utility::GSL_error("Vector::resize failed to allocate memory");
208    all(init_value);
209  }
210
211
212  void swap(Vector& v, Vector& w)
213  {
214    assert(v.gsl_vector_p()); assert(w.gsl_vector_p());
215    int status=gsl_vector_swap(v.gsl_vector_p(),w.gsl_vector_p());
216    if (status)
217      throw utility::GSL_error(std::string("swap(Vector&,Vector&)",status));
218  }
219
220
221  const Vector& Vector::operator=( const VectorBase& other )
222  {
223    return assign(other);
224  } 
225
226
227  const Vector& Vector::operator=( const Vector& other )
228  {
229    return assign(other);
230  }
231
232
233}}} // of namespace utility, yat, and thep
Note: See TracBrowser for help on using the repository browser.