vector3.h

#ifndef __VECTOR3_H
#define __VECTOR3_H

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

template <class REAL> struct vector3{
public:
	REAL x, y, z;
	vector3(){
		x = y = z = REAL(0);
	}
	vector3(const REAL &r){
		x = y = z = r;
	}
	vector3(const REAL &_x, const REAL &_y, const REAL &_z){
		x = _x;  y = _y;  z = _z;
	}
	vector3(const REAL *p){
		x = p[0]; y = p[1]; z = p[2];
	}
	~vector3(){}

	REAL &operator [](int i){
	  return (&x)[i];
	}
	template <class real> 
		operator vector3<real> () const {
			return vector3<real> (real(x), real(y), real(z));
		}
	operator REAL *(){
		return &x;
	}
	operator const REAL *() const{
		return &x;
	}
	REAL (*toPointer())[3]{
		return (REAL (*)[3])&x;
	}
	const REAL (*toConstPointer())[3]{
		return (const REAL (*)[3])&x;
	}
	typedef REAL (*pArrayOfReal3)[3];
	typedef const REAL (*cpArrayOfReal3)[3];
	operator pArrayOfReal3(){
		return toPointer();
	}
	operator cpArrayOfReal3(){
		return toConstPointer();
	}

	void outv(std::ostream &ofs = std::cout) const{
		ofs << "(" << x << ", " << y << ", " << z << ")" << std::endl;
	}
	bool are_numbers () const{
		// returns false if *this has (a) NaN member(s)
		return (norm2() >= REAL(0));
	}

	REAL norm2() const{
		return (*this)*(*this);
	}
	REAL abs() const{
		return std::sqrt(norm2());
	}

	friend std::ostream &operator << (std::ostream &ofs, const vector3<REAL> &v){
		ofs << v.x << " " << v.y << " " << v.z;
		return ofs;
	}
	friend std::istream &operator >> (std::istream &ifs, vector3<REAL> &v){
		ifs >> v.x >> v.y >> v.z;
		return ifs;
	}
	const vector3<REAL> operator + (const vector3<REAL> &v) const{
		return vector3<REAL> (x+v.x, y+v.y, z+v.z);
	}
	const vector3<REAL> operator - (const vector3<REAL> &v) const{
		return vector3<REAL> (x-v.x, y-v.y, z-v.z);
	}
	const vector3<REAL> operator * (const REAL &s) const{
		return vector3<REAL> (x*s, y*s, z*s);
	}
	friend const vector3<REAL> operator * (const REAL &s, const vector3<REAL> &v){
		return v*s;
	}
	// dot product
	const REAL operator * (const vector3<REAL> &v) const{
		return (x*v.x + y*v.y + z*v.z);
	}
	const vector3<REAL> cross(const vector3<REAL> &v) const 
	{
		return vector3<REAL>(y*v.z - z*v.y, z*v.x - x*v.z, x*v.y - y*v.x);
	}
	// vector product
	const vector3<REAL> operator % (const vector3<REAL> &v) const{
    return vector3<REAL> (
        y*v.z - z*v.y, 
        z*v.x - x*v.z, 
        x*v.y - y*v.x);
  }
  const vector3<REAL> operator / (const REAL &s) const{
    REAL r = REAL(1)/s;
    return (*this)*r;
  }
  const vector3<REAL> operator = (const vector3<REAL> &v){
    x = v.x; y=v.y; z=v.z;
    return *this;
  }

  const vector3<REAL> operator - () const {
    return vector3<REAL> (-x, -y, -z);
  }
  const vector3<REAL> &operator += (const vector3<REAL> &v){
    *this = *this + v;
    return *this;
  }
  const vector3<REAL> &operator -= (const vector3<REAL> &v){
    *this = *this - v;
    return *this;
  }
  const vector3<REAL> &operator *= (const REAL &s){
    *this = *this * s;
    return *this;
  }
  const vector3<REAL> &operator /= (const REAL &s){
    *this = *this / s;
    return *this;
  }

  friend const vector3<REAL> maxeach (const vector3<REAL> &a, const vector3<REAL> &b){
    return vector3<REAL> (std::max(a.x, b.x), std::max(a.y, b.y), std::max(a.z, b.z));
  }
  friend const vector3<REAL> mineach (const vector3<REAL> &a, const vector3<REAL> &b){
    return vector3<REAL> (std::min(a.x, b.x), std::min(a.y, b.y), std::min(a.z, b.z));
  }
  const vector3<REAL> abseach() const {
    return vector3<REAL> (std::fabs(x), std::fabs(y), std::fabs(z));
  }
};

typedef double real;
typedef vector3<real> vec3;

#endif //  __VECTOR3_H