Implementation of Nd vectors using valarray

This is my code which implements Nd vectors. I chose valarray as a base of the class. Is this a good idea or I simply need to make a typedef std::valarray gvector instead of writing a class? Please review this code.

//gvector.h
#ifndef GVECTOR_H_
#define GVECTOR_H_

#include <cmath>
#include <fstream>
#include <iostream>
#include <valarray>

double sqr(const double &); //square

typedef std::valarray < double > coordinates;

std::ifstream & operator >> (std::ifstream &, coordinates &);
std::ostream &  operator << (std::ostream &,  const coordinates &);

class gvector;
std::ifstream & operator >> (std::ifstream &, gvector &);
std::ostream &  operator << (std::ostream &,  const gvector &);

class gvector {

friend std::ifstream & operator >> (std::ifstream &, gvector &);
friend std::ostream &  operator << (std::ostream &,  const gvector &);

coordinates crds;

public:

gvector(): crds() {};

gvector(const size_t &      d): crds(d)          {}; //set dimension
gvector(const coordinates & c): crds(c)          {};
gvector(const gvector &     v): crds(v.coords()) {};

//Initialize gvector by its start-point and end-point:
gvector(
const coordinates & start,
const coordinates & end
): crds(end - start) {};

size_t      dim()    const; //show dimension
coordinates coords() const; //show
double      len()    const; //show length

void set_dim(const size_t &);
void set_coords(const coordinates &);
void set_coords(const coordinates &, const coordinates &); //Change start-point and end-point of gvector

double & operator [] (const size_t &);
double const & operator [] (const size_t &) const;

//Opposite vector
gvector operator - () const;

gvector & operator =  (const gvector &);
gvector & operator += (const gvector &);
gvector & operator -= (const gvector &);
gvector & operator *= (const double &);

};

inline size_t      gvector::dim()    const { return crds.size(); }
inline coordinates gvector::coords() const { return crds;        }

inline double gvector::len() const {
return sqrt(crds.apply(sqr).sum());
}

inline void gvector::set_dim(const size_t & d) {
if (crds.size() != d) crds.resize(d);
}

inline void gvector::set_coords(const coordinates & c) {
size_t csz = c.size();
if (crds.size() != csz) crds.resize(csz);
crds = c;
}

inline void gvector::set_coords(const coordinates & start, const coordinates & end) {
size_t  ssz = start.size(),
esz = end.size();
if (ssz == esz && ssz == crds.size()) crds = end - start;
}

inline double & gvector::operator [] (const size_t & idx) {
return crds[idx];
}

inline double const & gvector::operator [] (const size_t & idx) const {
return crds[idx];
}

inline gvector gvector::operator - () const { return gvector(-crds); }

inline gvector & gvector::operator = (const gvector & v) {
size_t vdim = v.dim();
if (crds.size() != vdim) crds.resize(vdim);
crds = v.coords();
return *this;
}

inline gvector & gvector::operator += (const gvector & v) {
if (crds.size() == v.dim()) crds += v.coords();
return *this;
}

inline gvector & gvector::operator -= (const gvector & v) {
if (crds.size() == v.dim()) crds -= v.coords();
return *this;
}

inline gvector & gvector::operator *= (const double & scalar) {
crds *= scalar;
return *this;
}

gvector operator + (gvector, const gvector &);
gvector operator - (gvector, const gvector &);
gvector operator * (gvector, const double &);
gvector operator * (const double &, gvector);

#endif /* GVECTOR_H_ */


//gvector.cpp

#include "gvector.h"
#include <iostream>
#include <fstream>

double sqr(const double & x) { return x * x; }

std::ostream & operator << (std::ostream & os, const coordinates & c) {
if (c.size()) {
os << '(' << c[0];
for (size_t i = 1; i < c.size(); ++i) {
os << ", " << c[i];
}
os << ')';
}
return os;
}

std::ifstream & operator >> (std::ifstream & ifs, coordinates & c) {
for (size_t i = 0; i < c.size(); ++i) {
ifs >> c[i];
}
return ifs;
}

std::ostream & operator << (std::ostream & os, const gvector & v) {
os << v.crds;
return os;
}

std::ifstream & operator >> (std::ifstream & ifs, gvector & v) {
ifs >> v.crds;
return ifs;
}

gvector operator + (gvector lft, const gvector & rht) {
return lft += rht;
}

gvector operator - (gvector lft, const gvector & rht) {
return lft -= rht;
}

gvector operator * (gvector v, const double & scalar) {
return v *= scalar;
}

gvector operator * (const double & scalar, gvector v) {
return v *= scalar;
}


std::valarray by itself only acts like a 1D array. The <valarray> header includes a few ancillary classes to do this though. More importantly, valarray itself includes direct support for those classes.

For your purposes, the most important of these are probably std::slice and std::gslice. std::slice allows you to address a slice of the array and address it (about) like you'd index along a dimension of a multi-dimensional array.

For example, let's assume we have an array of 28 items that we want to treat as a 2D array (with dimensions 4x7). While it's probably not entirely obvious from the standard1 std::slice makes it fairly easy to do this. Let's say we want to walk across it, and look at one "column" of the 2D array at a time. So, we're going to create 7 columns (one at a time) each containing 4 items.

Now, the trick to this: the slices you create this way act as objects in themselves. You don't get access to the individual elements in the slice; you get access to the slice as a whole. To assign one value to all the elements in the slice, you do not walk through the slice with a for loop and assign to each element individually. Rather, you assign to the entire slice at once (which will assign that value to all the elements in the slice).

If you want to assign separate/unrelated values to the items in the slice, you can do that--but you need to use some other slice or valarray as the source. For example, let's create a valarray of 28 items, which we'll treat as a 2D array of 4x7 (4 rows of 7 columns each). For the sake of argument, we'll start with a valarray of 4 items holding 1, 2, 3, 4. We'll then assign that to the columns in the matrix. To keep things from getting too boring, we'll multiply that by the column number before assigning it, so the first row should count across by 1's, the second row by 2's, the third row by 3's and the fourth row by 4's.

std::valarray<int> data(rows * cols);
std::valarray<int> filler = { 1, 2, 3, 4 };

for (int i = 0; i < cols; i++)
data[std::slice(i, rows, cols)] = filler * i;

// show the data:
for (int i = 1; i <= t.size(); i++) {
std::cout << t[i-1] << "\t";
if (i % cols == 0)
std::cout << "\n";
}


The result should look like this:

0       1       2       3       4       5       6
0       2       4       6       8       10      12
0       3       6       9       12      15      18
0       4       8       12      16      20      24


As we can see, this fits the expected results, indicating that our slice has allowed us to address columns of the array, exactly as intended.

So, the bottom line: no, you can't (quite) just use std::valarray to do N-dimensional array access--but you can do it with valarray + the other classes included in it header. As such, I'm not sure your gvector class really adds any capability that isn't already present.

1. Okay, almost nothing in the C++ standard is really obvious, but this is probably even harder to figure out than most things.