nemo/Vector.h
Go to the documentation of this file.00001 #ifndef _NEMO_VECTOR_H_ 00002 #define _NEMO_VECTOR_H_ 00003 00004 00005 #include "IntrusiveCowSupport.h" 00006 #include "IntrusiveCounted.h" 00007 #include "Collection.h" 00008 #include "detail/nemo-cpp11.h" 00009 #include "detail/nemo-errors.h" 00010 #include "Random.h" 00011 00012 00013 00014 namespace nemo { 00015 00016 00017 class VectorDimension; 00018 00019 static VectorDimension dim(unsigned int d); 00020 00025 class VectorDimension { 00026 public: 00027 inline const unsigned int& get() const {return d;} 00028 friend VectorDimension dim(unsigned int d); 00029 operator unsigned int () const { 00030 return d; 00031 } 00032 private: 00033 VectorDimension(unsigned int dd) : d(dd) {} 00034 //VectorDimension(VectorDimension &o) : d(o.d) {} 00035 //VectorDimension& operator=(const VectorDimension& o){d = o.d; return *this;} 00036 unsigned int d; 00037 }; 00038 00042 static inline VectorDimension dim(unsigned int d){ 00043 return VectorDimension(d); 00044 } 00045 00046 template <typename T> 00047 class VectorTranspose; 00048 00067 template <class T> 00068 class MathVector : public Collection<T> { 00069 00070 00071 #ifndef NEMO_PARSED_BY_DOXYGEN 00072 00073 #define _NEMO_VECTOR_CHECK_INDEX(index) \ 00074 if(data->dim <= index){ \ 00075 _NEMO_ARGUMENT_ERROR("MathVector: index out of range") \ 00076 } 00077 00078 private: 00079 friend class VectorTranspose<T>; 00080 00081 public: 00082 //class CWise; 00083 friend class CWise; // CWise needs to access CWiseOperation classes 00084 #endif 00085 00086 private: 00087 00088 // ============================================================================= 00089 // Data management 00090 // ============================================================================= 00091 class VectorData; 00092 00093 #ifndef NEMO_PARSED_BY_DOXYGEN 00094 00098 class VectorDataDeleter { 00099 public: 00100 VectorDataDeleter(const void *ptr){ 00101 // do NOT just delete the C++ object! 00102 // object plus data array must be deleted alltogether by means of 'free', 00103 // accordingly to 'malloc'-creation inside VectorData::create() 00104 free((VectorData*)ptr); 00105 } 00106 }; 00107 00115 class VectorData : public IntrusiveCounted<VectorData,VectorDataDeleter> { 00116 public: 00118 // byte 0-3: reference counter, inherited from IntrusiveCounted 00119 // byte 4-7: dimension count of the vector: 00120 boost::uint32_t dim; 00121 // byte 8-11: pointer to data array (on 64 bit it's byte 8-15!) 00122 T *values; 00123 // byte 12-X: data array. not part of the C++ object, but allocated 00124 // right here by the create()-method (on 64 bit it's byte 16-X!) 00125 00127 static VectorData* create(const unsigned int dimension){ 00128 // alloc VectorData structure, and directly behind it the space for the data-array 00129 VectorData* ptr = (VectorData*)malloc( sizeof(VectorData) + dimension*sizeof(T) ); 00130 ptr->initRefCount(); 00131 ptr->dim = dimension; 00132 // the values now lie at ( &values (as void for single byte ptr arithm.) + sizeof(T*)bytes ) 00133 ptr->values = (T*)((boost::uint8_t*)&(ptr->values)+sizeof(T*)); 00134 return ptr; 00135 } 00136 private: 00137 // these functions MUST NOT be used! use VectorData::create/CustomCopyAllocator/VectorDataDeleter instead 00138 VectorData(){_NEMO_LOGIC_ERROR("MathVector::VectorData: ctor must not be called");} 00139 VectorData(const VectorData&){_NEMO_LOGIC_ERROR("MathVector::VectorData: copy-ctor must not be called");} 00140 ~VectorData(){_NEMO_LOGIC_ERROR("MathVector::VectorData: dtor must not be called");} 00141 }; 00142 00147 class CustomCopyAllocator { 00148 public: 00149 static VectorData* copyAllocate(const VectorData* vd) { 00150 // the normal copy-constructor MUST NOT be used here, 00151 // since it would only re-create the object, but not the data array behind it! 00152 VectorData* ptr = VectorData::create(vd->dim); 00153 for(unsigned int i=0; i<vd->dim; i++){ 00154 ptr->values[i] = vd->values[i]; 00155 } 00156 return ptr; 00157 } 00158 }; 00159 00160 typedef IntrusiveCowSupport<VectorData,CustomCopyAllocator> VectorDataPtr; 00161 #endif 00162 00163 public: 00164 // ============================================================================= 00165 // Construction 00166 // ============================================================================= 00167 00171 explicit MathVector() 00172 : data(VectorData::create(0)) 00173 {} 00174 00182 explicit MathVector(const VectorDimension &dim) 00183 : data(VectorData::create(dim.get())) 00184 { 00185 for(unsigned int i=0; i<dim.get(); i++){ 00186 data->values[i] = T(); 00187 } 00188 } 00189 00197 explicit MathVector(const VectorDimension &dim, const T &val) 00198 : data(VectorData::create(dim.get())) 00199 { 00200 for(unsigned int i=0; i<dim.get(); i++){ 00201 data->values[i] = val; 00202 } 00203 } 00204 00205 // FIXME maybe the template should be removed by typesafe Generator<T>&, avoid possible confusion with Vec(dim,val) 00216 template <typename generatorType> 00217 explicit MathVector(const VectorDimension &dim, generatorType gen) 00218 : data(VectorData::create(dim.get())) 00219 { 00220 for(unsigned int i=0; i<dim.get(); i++){ 00221 data->values[i] = gen(); 00222 } 00223 } 00224 00228 explicit MathVector(const T &x) 00229 : data(VectorData::create(1)) 00230 { 00231 data->values[0] = x; 00232 } 00236 explicit MathVector(const T &x, const T &y) 00237 : data(VectorData::create(2)) 00238 { 00239 data->values[0] = x; 00240 data->values[1] = y; 00241 } 00245 explicit MathVector(const T &x, const T &y, const T &z) 00246 : data(VectorData::create(3)) 00247 { 00248 data->values[0] = x; 00249 data->values[1] = y; 00250 data->values[2] = z; 00251 } 00252 00253 #ifdef _NEMO_HAS_INIT_LISTS_ 00254 MathVector(const std::initializer_list<T> &init) 00255 : data(VectorData::create(init.size())) 00256 { 00257 unsigned int i=0; 00258 for(const T *p = init.begin (); p != init.end (); ++p){ 00259 data->values[i++] = *p; 00260 } 00261 } 00262 #endif 00263 00264 MathVector<T>(const MathVector<T> &v) 00265 : data(v.data) {} 00266 00267 MathVector<T>& operator=(const MathVector<T> &v){ 00268 this->data = v.data; 00269 return *this; 00270 } 00271 00278 template <typename D> 00279 static MathVector<T> fromArray(const VectorDimension &dim, const D* vals) { 00280 MathVector<T> vec(dim); 00281 for(unsigned int i=0; i<dim.get(); i++){ 00282 vec.data->values[i] = vals[i]; 00283 } 00284 return vec; 00285 } 00286 00287 // ============================================================================= 00288 // Access 00289 // ============================================================================= 00290 00294 //unsigned int dimension() const { 00295 // return data->dim; 00296 //} 00297 inline VectorDimension dimension() const { 00298 return dim(data->dim); 00299 } 00300 00305 const T& operator[](const unsigned int index) const{ 00306 _NEMO_VECTOR_CHECK_INDEX(index); 00307 return data->values[index]; 00308 } 00309 00310 #ifndef NEMO_PARSED_BY_DOXYGEN 00311 class ReadWriteRef { 00312 friend class MathVector; 00313 friend class MathVector::Iterator; 00314 private: 00315 ReadWriteRef(VectorDataPtr &d, const unsigned int i) : data(d), index(i) {} 00316 public: 00317 ReadWriteRef& operator=(const T &value){ 00318 // write access 00319 data.detach(); 00320 data->values[index] = value; 00321 return *this; 00322 } 00323 ReadWriteRef& operator=(ReadWriteRef &value){ 00324 // write access 00325 data.detach(); 00326 data->values[index] = (T)value; 00327 return *this; 00328 } 00329 00330 operator const T() const{return data->values[index];} 00331 00332 #define _nemo_mathvector_ref_define_compound_(OP) \ 00333 ReadWriteRef& operator OP (const T& value){\ 00334 data.detach();\ 00335 data->values[index] OP value;\ 00336 return *this;\ 00337 } 00338 _nemo_mathvector_ref_define_compound_(+=) 00339 _nemo_mathvector_ref_define_compound_(-=) 00340 _nemo_mathvector_ref_define_compound_(*=) 00341 _nemo_mathvector_ref_define_compound_(/=) 00342 _nemo_mathvector_ref_define_compound_(%=) 00343 _nemo_mathvector_ref_define_compound_(&=) 00344 _nemo_mathvector_ref_define_compound_(|=) 00345 _nemo_mathvector_ref_define_compound_(^=) 00346 _nemo_mathvector_ref_define_compound_(<<=) 00347 _nemo_mathvector_ref_define_compound_(>>=) 00348 00349 #define _nemo_mathvector_ref_define_compare_(OP) \ 00350 bool operator OP (const typename MathVector<T>::ReadWriteRef& rhs) const{\ 00351 return data->values[index] OP rhs.data->values[rhs.index];\ 00352 } 00353 _nemo_mathvector_ref_define_compare_(==) 00354 _nemo_mathvector_ref_define_compare_(!=) 00355 _nemo_mathvector_ref_define_compare_(<=) 00356 _nemo_mathvector_ref_define_compare_(<) 00357 _nemo_mathvector_ref_define_compare_(>=) 00358 _nemo_mathvector_ref_define_compare_(>) 00359 00360 /* gives wearnings not to return ReadWriteRef& 00361 * pragma only works globally *argh* 00362 *#pragma GCC diagnostic ignored "-Weffc++" 00363 * standard signature would be 00364 FOO& FOO::operator ++(); 00365 FOO FOO::operator ++(int); 00366 * but prefix operator must not return ref to itself, 00367 * that would behave like the postfix op. Rather, the type 00368 * shuold an object of type T, but not writable (don't return T&!!!) 00369 */ 00370 /* 00371 // prefix 00372 T operator++(){ 00373 data.detach(); 00374 return ++(data->values[index]); 00375 } 00376 // postfix 00377 T operator++(int){ 00378 data.detach(); 00379 return (data->values[index])++; 00380 } 00381 00382 // prefix 00383 T operator--(){ 00384 data.detach(); 00385 return --(data->values[index]); 00386 } 00387 // postfix 00388 T operator--(int){ 00389 data.detach(); 00390 return (data->values[index])--; 00391 } 00392 */ 00393 00394 // this is a 'specialization' of std::swap that can be found by argument-dependent lookup. 00395 // this one is necessary for STL write operations to work correctly, since the standard swap 00396 // implementation does not work correctly on ReadWriteRef, since it already has reference semantics 00397 // note that arguments are passed by value (copy-construction does no internal write), 00398 // in order to permit passing of temp. values from operator[] 00399 friend inline void swap( typename nemo::MathVector<T>::ReadWriteRef a, typename nemo::MathVector<T>::ReadWriteRef b ){ 00400 a.data.detach(); 00401 b.data.detach(); 00402 using std::swap; 00403 swap( a.data->values[a.index], b.data->values[b.index] ); 00404 } 00405 friend inline T max( const typename nemo::MathVector<T>::ReadWriteRef& a, const typename nemo::MathVector<T>::ReadWriteRef& b ){ 00406 using std::max; 00407 return max( a.data->values[a.index], b.data->values[b.index] ); 00408 } 00409 friend inline T max( const typename nemo::MathVector<T>::ReadWriteRef& a, T b ){ 00410 using std::max; 00411 return max( a.data->values[a.index], b ); 00412 } 00413 friend inline T max( T a, const typename nemo::MathVector<T>::ReadWriteRef& b ){ 00414 using std::max; 00415 return max( a, b.data->values[b.index] ); 00416 } 00417 friend inline T min( const typename nemo::MathVector<T>::ReadWriteRef& a, const typename nemo::MathVector<T>::ReadWriteRef& b ){ 00418 using std::min; 00419 return min( a.data->values[a.index], b.data->values[b.index] ); 00420 } 00421 friend inline T min( const typename nemo::MathVector<T>::ReadWriteRef& a, T b ){ 00422 using std::min; 00423 return min( a.data->values[a.index], b ); 00424 } 00425 friend inline T min( T a, const typename nemo::MathVector<T>::ReadWriteRef& b ){ 00426 using std::min; 00427 return min( a, b.data->values[b.index] ); 00428 } 00429 00430 private: 00431 VectorDataPtr &data; 00432 unsigned int index; 00433 }; 00434 #endif 00435 00436 00442 inline ReadWriteRef operator[](const unsigned int index){ 00443 _NEMO_VECTOR_CHECK_INDEX(index); 00444 return ReadWriteRef(data, index); // FIXME assignment/copy??? 00445 } 00446 00447 00454 const MathVector<T> subVector(const unsigned int fromIndex, const VectorDimension &dim) const { 00455 if(data->dim == 0){ 00456 if(fromIndex!=0 || dim.get()!=0){ 00457 _NEMO_ARGUMENT_ERROR( 00458 "zero-dimensional vector can only be devided with subVector(0,0)"); 00459 } 00460 return *this; 00461 } 00462 00463 _NEMO_VECTOR_CHECK_INDEX(fromIndex); 00464 if(0 < dim.get()){ 00465 _NEMO_VECTOR_CHECK_INDEX(fromIndex+dim.get()-1); 00466 } 00467 00468 MathVector<T> sub(dim); 00469 for(unsigned int i=0; i<dim.get(); i++){ 00470 sub.data->values[i] = data->values[i+fromIndex]; 00471 } 00472 return sub; 00473 } 00474 00475 void setSubVector(const unsigned int fromIndex, const MathVector<T> &vec){ 00476 _NEMO_VECTOR_CHECK_INDEX(fromIndex+vec.dimension()-1); 00477 00478 data.detach(); 00479 for(unsigned int i=0; i<vec.dimension(); i++){ 00480 data->values[fromIndex+i] = vec.data->values[i]; 00481 } 00482 } 00483 00484 // ============================================================================= 00485 // Comparison 00486 // ============================================================================= 00491 bool operator==(const MathVector<T> &r) const{ 00492 if(data->dim != r.data->dim){ 00493 return false; 00494 } 00495 for(unsigned int i=0; i<data->dim; i++){ 00496 if(data->values[i] != r.data->values[i]){ 00497 return false; 00498 } 00499 } 00500 // else 00501 return true; 00502 } 00503 00507 bool operator!=(const MathVector<T> &r) const{ 00508 return !operator==(r); 00509 } 00510 00511 // ============================================================================= 00512 // Collection API 00513 // ============================================================================= 00514 00515 unsigned int numElements() const { 00516 return dimension(); 00517 } 00518 T getElement(unsigned int elementIndex) const{ 00519 return (*this)[elementIndex]; 00520 } 00521 void setElement(unsigned int elementIndex, T value){ 00522 (*this)[elementIndex] = value; 00523 } 00524 00525 void transform(Mapping<T,T> m){ 00526 for(unsigned int i=0; i<dimension(); i++){ 00527 (*this)[i] = m((*this)[i]); 00528 } 00529 } 00530 00531 template <typename OutType> 00532 MathVector<OutType> map(Mapping<T,OutType> m) const { 00533 MathVector<OutType> out(dim(this->dimension())); 00534 for(unsigned int i=0; i<dimension(); i++){ 00535 out[i] = m((*this)[i]); 00536 } 00537 return out; 00538 } 00539 00540 // ============================================================================= 00541 // STL compliant iterator interface 00542 // ============================================================================= 00543 00544 #ifndef NEMO_PARSED_BY_DOXYGEN 00545 class Iterator{ 00546 private: 00547 Iterator(const MathVector<T> *vecP, int indexP) : vec(vecP), index(indexP), rwRef(const_cast<MathVector<T>*>(vecP)->data,(unsigned int)indexP) {} 00548 public: 00549 Iterator(const Iterator& it) : vec(it.vec), index(it.index), rwRef(const_cast<MathVector<T>*>(it.vec)->data,(unsigned int)it.index) {} 00550 Iterator& operator=(const Iterator &it){ 00551 if(this==&it){return *this;} 00552 this->vec = it.vec; 00553 this->index = it.index; 00554 this->rwRef.data = it.rwRef.data; 00555 this->rwRef.index = it.rwRef.index; 00556 return *this; 00557 } 00558 00559 typedef int difference_type; 00560 typedef T value_type; 00561 typedef ReadWriteRef* pointer; 00562 typedef ReadWriteRef reference; 00563 //typedef std::bidirectional_iterator_tag iterator_category; 00564 typedef std::random_access_iterator_tag iterator_category; 00565 00566 ReadWriteRef& operator*(){ 00567 if(index<0 || (int)vec->dimension()<=index){_NEMO_LOGIC_ERROR("Iterator out of range");} 00568 rwRef.index = (unsigned int)index; 00569 return rwRef; 00570 //return (*const_cast<MathVector<T>*>(vec))[(unsigned int)index]; 00571 } 00572 T operator*() const{ 00573 if(index<0 || (int)vec->dimension()<=index){_NEMO_LOGIC_ERROR("Iterator out of range");} 00574 return vec->getElement((unsigned int)index); 00575 } 00576 ReadWriteRef& operator[](difference_type offset){ 00577 if(index+offset<0 || (int)vec->dimension()<=index+offset){_NEMO_LOGIC_ERROR("Iterator out of range");} 00578 rwRef.index = (unsigned int)index+offset; 00579 return rwRef; 00580 //return (*const_cast<MathVector<T>*>(vec))[(unsigned int)index+offset]; 00581 } 00582 T operator[](difference_type offset) const{ 00583 if(index+offset<0 || (int)vec->dimension()<=index+offset){_NEMO_LOGIC_ERROR("Iterator out of range");} 00584 return vec->getElement((unsigned int)index+offset); 00585 } 00586 Iterator& operator++() const{ 00587 index++; 00588 return const_cast<Iterator&>(*this); 00589 } 00590 Iterator operator++(int) const{ 00591 Iterator copy = *this; 00592 index++; 00593 return copy; 00594 } 00595 Iterator& operator--() const{ 00596 index--; 00597 return const_cast<Iterator&>(*this); 00598 } 00599 Iterator operator--(int) const{ 00600 Iterator copy = *this; 00601 index--; 00602 return copy; 00603 } 00604 Iterator operator+(difference_type off) const { 00605 Iterator it = *this; 00606 it.index += off; 00607 return it; 00608 } 00609 Iterator operator-(difference_type off) const { 00610 Iterator it = *this; 00611 it.index -= off; 00612 return it; 00613 } 00614 difference_type operator-(const Iterator& rhs) const { 00615 return index - rhs.index; 00616 } 00617 bool operator==(const Iterator &it) const{ 00618 return (vec==it.vec) && (index==it.index); 00619 } 00620 bool operator!=(const Iterator &it) const{ 00621 return !(*this==it); 00622 } 00623 bool operator<=(const Iterator &it) const{ 00624 return index<=it.index; 00625 } 00626 bool operator<(const Iterator &it) const{ 00627 return index<it.index; 00628 } 00629 bool operator>=(const Iterator &it) const{ 00630 return index>=it.index; 00631 } 00632 bool operator>(const Iterator &it) const{ 00633 return index>it.index; 00634 } 00635 private: 00636 friend class MathVector<T>; 00637 const MathVector<T> *vec; 00638 mutable int index; 00639 mutable ReadWriteRef rwRef; 00640 }; 00641 #endif 00642 00646 Iterator begin(){ 00647 return Iterator(this,0); 00648 } 00652 Iterator end(){ 00653 return Iterator(this,(int)data->dim); 00654 } 00655 00656 const Iterator begin() const{ 00657 return Iterator(this,0); 00658 } 00659 const Iterator end() const{ 00660 return Iterator(this,(int)data->dim); 00661 } 00662 00663 // ============================================================================= 00664 // Vector-Vector operations 00665 // ============================================================================= 00666 00667 private: 00668 #ifndef NEMO_PARSED_BY_DOXYGEN 00669 class VectorAdd{ 00670 public: 00671 static void evaluate(const MathVector<T> &v1, const MathVector<T> &v2, MathVector<T> &result){ 00672 for(unsigned int i=0; i<result.data->dim; i++){ 00673 result.data->values[i] = v1.data->values[i] + v2.data->values[i]; 00674 } 00675 } 00676 }; 00677 class VectorDiff{ 00678 public: 00679 static void evaluate(const MathVector<T> &v1, const MathVector<T> &v2, MathVector<T> &result){ 00680 for(unsigned int i=0; i<result.data->dim; i++){ 00681 result.data->values[i] = v1.data->values[i] - v2.data->values[i]; 00682 } 00683 } 00684 }; 00685 class VectorMultScalar{ 00686 public: 00687 static void evaluate(const MathVector<T> &v, const T &scalar, MathVector<T> &result){ 00688 for(unsigned int i=0; i<result.data->dim; i++){ 00689 result.data->values[i] = v.data->values[i] * scalar; 00690 } 00691 } 00692 }; 00693 class VectorCWiseMult{ 00694 public: 00695 static void evaluate(const MathVector<T> &v1, const MathVector<T> &v2, MathVector<T> &result){ 00696 for(unsigned int i=0; i<result.data->dim; i++){ 00697 result.data->values[i] = v1.data->values[i] * v2.data->values[i]; 00698 } 00699 } 00700 }; 00701 class VectorCWiseDiv{ 00702 public: 00703 static void evaluate(const MathVector<T> &v1, const MathVector<T> &v2, MathVector<T> &result){ 00704 for(unsigned int i=0; i<result.data->dim; i++){ 00705 result.data->values[i] = v1.data->values[i] / v2.data->values[i]; 00706 } 00707 } 00708 }; 00709 class VectorCWiseMax{ 00710 public: 00711 static void evaluate(const MathVector<T> &v1, const MathVector<T> &v2, MathVector<T> &result){ 00712 for(unsigned int i=0; i<result.data->dim; i++){ 00713 result.data->values[i] = 00714 (v1.data->values[i] < v2.data->values[i]) ? 00715 v2.data->values[i] : v1.data->values[i]; 00716 } 00717 } 00718 }; 00719 class VectorCWiseMin{ 00720 public: 00721 static void evaluate(const MathVector<T> &v1, const MathVector<T> &v2, MathVector<T> &result){ 00722 for(unsigned int i=0; i<result.data->dim; i++){ 00723 result.data->values[i] = 00724 (v1.data->values[i] < v2.data->values[i]) ? 00725 v1.data->values[i] : v2.data->values[i]; 00726 } 00727 } 00728 }; 00729 #endif 00730 00731 /*template <class Op> 00732 class VectorBinaryResult { 00733 public: 00734 VectorBinaryResult(const MathVector<T> &v1P, const MathVector<T> &v2P) 00735 : v1(v1P), v2(v2P) 00736 {} 00737 00738 operator MathVector<T>() const { 00739 MathVector<T> r(dim(v1.dimension())); 00740 Op::evaluate(v1, v2, r); 00741 return r; 00742 } 00743 00744 private: 00745 const MathVector<T> &v1; 00746 const MathVector<T> &v2; 00747 }; 00748 00749 template <typename Op> 00750 MathVector(const VectorBinaryResult<Op> &result) 00751 : data(new VectorData(result.v1.data->dim)) 00752 { 00753 Op::evaluate(result.v1, result.v2, *this); 00754 } 00755 template <typename Op> 00756 MathVector<T>& operator=(const VectorBinaryResult<Op> &result){ 00757 if(data->dim == result.v1.data->dim){ 00758 // TODO this should be uninitialized? 00759 data.detach(); 00760 } 00761 else{ 00762 data = new VectorData(result.v1.data->dim); 00763 } 00764 assert(data->dim == result.v1.data->dim); 00765 00766 Op::evaluate(result.v1, result.v2, *this); 00767 00768 return *this; 00769 }*/ 00770 00771 00772 public: 00776 /*const VectorBinaryResult<VectorAdd> operator+(MathVector<T> &r){ 00777 checkDimensionMatch(r); 00778 return VectorBinaryResult<VectorAdd>(*this, r); 00779 }*/ 00780 MathVector<T> operator+(const MathVector<T> &r) const{ 00781 checkDimensionMatch(r); 00782 MathVector<T> result(dim(data->dim)); 00783 VectorAdd::evaluate(*this, r, result); 00784 return result; 00785 } 00789 MathVector<T>& operator+=(const MathVector<T> &r){ 00790 checkDimensionMatch(r); 00791 data.detach(); 00792 VectorAdd::evaluate(*this, r, *this); 00793 return *this; 00794 } 00795 00799 /*const VectorBinaryResult<VectorDiff> operator-(MathVector<T> &r) const{ 00800 checkDimensionMatch(r); 00801 return VectorBinaryResult<VectorDiff>(*this, r); 00802 }*/ 00803 MathVector<T> operator-(const MathVector<T> &r) const{ 00804 checkDimensionMatch(r); 00805 MathVector<T> result(dim(data->dim)); 00806 VectorDiff::evaluate(*this, r, result); 00807 return result; 00808 } 00812 MathVector<T>& operator-=(const MathVector<T> &r){ 00813 checkDimensionMatch(r); 00814 data.detach(); 00815 VectorDiff::evaluate(*this, r, *this); 00816 return *this; 00817 } 00818 00819 00823 friend MathVector<T> operator*(const MathVector<T> &vec, const T &scalar){ 00824 MathVector<T> result(dim(vec.dimension())); 00825 VectorMultScalar::evaluate(vec, scalar, result); 00826 return result; 00827 } 00831 friend MathVector<T> operator*(const T &scalar, const MathVector<T> &vec){ 00832 MathVector<T> result(dim(vec.dimension())); 00833 VectorMultScalar::evaluate(vec, scalar, result); 00834 return result; 00835 } 00839 MathVector<T>& operator*=(const T &scalar){ 00840 data.detach(); 00841 VectorMultScalar::evaluate(*this, scalar, *this); 00842 return *this; 00843 } 00844 00845 00851 MathVector<T> attach(const MathVector<T> &vec) const { 00852 //VectorData *newData = new VectorData(data->dim + vec.data->dim); 00853 MathVector<T> newVec(dim(data->dim + vec.data->dim)); 00854 for(unsigned int i=0; i<data->dim; i++){ 00855 newVec.data->values[i] = data->values[i]; 00856 } 00857 for(unsigned int i=0; i<vec.data->dim; i++){ 00858 newVec.data->values[data->dim+i] = vec.data->values[i]; 00859 } 00860 return newVec; 00861 } 00862 00863 // ============================================================================= 00864 // Vector unary operations 00865 // ============================================================================= 00866 private: 00867 #ifndef NEMO_PARSED_BY_DOXYGEN 00868 class VectorMinus{ 00869 public: 00870 static void evaluate(const MathVector<T> &v, MathVector<T> &result){ 00871 for(unsigned int i=0; i<result.data->dim; i++){ 00872 result.data->values[i] = - v.data->values[i]; 00873 } 00874 } 00875 }; 00876 #endif 00877 00878 public: 00879 MathVector<T> operator-() const{ 00880 MathVector<T> result(dim(data->dim)); 00881 VectorMinus::evaluate(*this, result); 00882 return result; 00883 } 00884 MathVector<T> operator+() const{ 00885 return *this; 00886 } 00887 00888 00889 // ============================================================================= 00890 // Vector-Vector component-wise operations 00891 // ============================================================================= 00892 00897 class CWise { 00898 friend class MathVector; 00899 private: 00900 CWise(const MathVector<T> &v) : vec(v) {} 00901 public: 00902 CWise() : vec() {} // needed for mapping operators! 00903 private: 00904 MathVector<T> vec; 00905 public: 00908 MathVector<T> operator*(const MathVector<T> &o) const { 00909 vec.checkDimensionMatch(o); 00910 MathVector<T> result(dim(vec.data->dim)); 00911 VectorCWiseMult::evaluate(vec, o, result); 00912 return result; 00913 } 00914 friend MathVector<T> operator*(const MathVector<T> &lhs, const CWise &rhs) { 00915 if(lhs.dimension() != rhs.vec.dimension()){ 00916 _NEMO_ARGUMENT_ERROR("MathVector::CWise::operator*: Dimension mismatch"); 00917 } 00918 MathVector<T> result(dim(lhs.dimension())); 00919 00920 //VectorCWiseMult::evaluate(lhs, rhs.vec, result); 00921 // FIXME the above call is not permitted by clang++ (private member, despite friend decl.) 00922 for(unsigned int i=0; i<result.dimension(); i++){ 00923 result[i] = lhs[i] * rhs.vec[i]; 00924 } 00925 00926 return result; 00927 } 00930 MathVector<T> operator/(const MathVector<T> &o) const { 00931 vec.checkDimensionMatch(o); 00932 MathVector<T> result(dim(vec.data->dim)); 00933 VectorCWiseDiv::evaluate(vec, o, result); 00934 return result; 00935 } 00936 friend MathVector<T> operator/(const MathVector<T> &lhs, const CWise &rhs) { 00937 if(lhs.dimension() != rhs.vec.dimension()){ 00938 _NEMO_ARGUMENT_ERROR("MathVector::CWise::operator/: Dimension mismatch"); 00939 } 00940 MathVector<T> result(dim(lhs.dimension())); 00941 00942 //VectorCWiseDiv::evaluate(lhs, rhs.vec, result); 00943 // FIXME the above call is not permitted by clang++ (private member, despite friend decl.) 00944 for(unsigned int i=0; i<result.dimension(); i++){ 00945 result[i] = lhs[i] / rhs.vec[i]; 00946 } 00947 return result; 00948 } 00951 MathVector<T> max(const MathVector<T> &o) const { 00952 vec.checkDimensionMatch(o); 00953 MathVector<T> result(dim(vec.data->dim)); 00954 VectorCWiseMax::evaluate(vec, o, result); 00955 return result; 00956 } 00959 MathVector<T> min(const MathVector<T> &o) const { 00960 vec.checkDimensionMatch(o); 00961 MathVector<T> result(dim(vec.data->dim)); 00962 VectorCWiseMin::evaluate(vec, o, result); 00963 return result; 00964 } 00965 operator const MathVector<T>() const {return vec;} 00966 }; 00967 00971 const CWise cwise() const { 00972 return CWise(*this); 00973 } 00974 00977 MathVector<T>& operator*=(const typename MathVector<T>::CWise &o){ 00978 checkDimensionMatch(o.vec); 00979 data.detach(); 00980 VectorCWiseMult::evaluate(*this, o.vec, *this); 00981 return *this; 00982 } 00985 MathVector<T>& operator/=(const typename MathVector<T>::CWise &o){ 00986 checkDimensionMatch(o.vec); 00987 data.detach(); 00988 VectorCWiseDiv::evaluate(*this, o.vec, *this); 00989 return *this; 00990 } 00991 00992 private: 00993 #ifndef NEMO_PARSED_BY_DOXYGEN 00994 template <typename OutType> 00995 class CWiseMapping : public MappingDefinition<MathVector<T>,MathVector<OutType> >{ 00996 public: 00997 MathVector<OutType> evaluate(MathVector<T> value) const{ 00998 return value.map(m); 00999 } 01000 CWiseMapping(Mapping<T,OutType> mP) : m(mP) {} 01001 private: 01002 Mapping<T,OutType> m; 01003 }; 01004 #endif 01005 01006 public: 01007 01012 template <typename OutType> 01013 static Mapping<MathVector<T>,MathVector<OutType> > cwise(Mapping<T,OutType> m){ 01014 return Mapping<MathVector<T>,MathVector<OutType> >(new CWiseMapping<OutType>(m)); 01015 } 01016 01017 // ============================================================================= 01018 // Vector-Vector operations with transpose logic 01019 // ============================================================================= 01020 01021 01025 VectorTranspose<T> transpose(){ 01026 return VectorTranspose<T>(*this); 01027 } 01028 01029 private: 01030 01031 template <class D> 01032 inline void checkDimensionMatch(const MathVector<D> &r) const{ 01033 if(data->dim != r.data->dim){ 01034 _NEMO_ARGUMENT_ERROR("MathVector: dimension mismatch"); 01035 } 01036 } 01037 01038 VectorDataPtr data; 01039 01040 }; 01041 01042 01046 template <typename T> 01047 class VectorTranspose { 01048 friend class MathVector<T>; 01049 private: 01050 #ifndef NEMO_PARSED_BY_DOXYGEN 01051 class VectorScalarProduct{ 01052 public: 01053 static void evaluate(const MathVector<T> &v1, const MathVector<T> &v2, T &result){ 01054 result = T(); 01055 for(unsigned int i=0; i<v1.data->dim; i++){ 01056 result += v1.data->values[i] * v2.data->values[i]; 01057 } 01058 } 01059 }; 01060 #endif 01061 01062 VectorTranspose(MathVector<T> &v) : vec(v) {} 01063 01064 MathVector<T> vec; 01065 public: 01066 VectorTranspose(const VectorTranspose &o) : vec(o.vec) {} 01067 VectorTranspose& operator=(const VectorTranspose &o){this->vec=o.vec; return *this;} 01068 01069 const T& operator[](const unsigned int index) const{ 01070 return vec[index]; 01071 } 01072 unsigned int dimension() const { 01073 return vec.dimension(); 01074 } 01078 T operator*(const MathVector<T> &o) const { 01079 vec.checkDimensionMatch(o); 01080 T result; 01081 VectorScalarProduct::evaluate(vec, o, result); 01082 return result; 01083 } 01084 01085 MathVector<T> transpose() const { 01086 return vec; 01087 } 01088 01089 bool operator==(const VectorTranspose<T> &r) const{ 01090 if(vec.data->dim != r.vec.data->dim){ 01091 return false; 01092 } 01093 for(unsigned int i=0; i<vec.data->dim; i++){ 01094 if(vec.data->values[i] != r.vec.data->values[i]){ 01095 return false; 01096 } 01097 } 01098 // else 01099 return true; 01100 } 01101 01105 bool operator!=(const VectorTranspose<T> &r) const{ 01106 return !operator==(r); 01107 } 01108 }; 01109 01110 01111 template <class T> 01112 std::ostream& operator<<(std::ostream& out, const MathVector<T> &vec){ 01113 01114 out << std::string("["); 01115 if(0 < vec.dimension()){ 01116 for(unsigned int i=0; i<vec.dimension()-1; i++){ 01117 out << vec[i] << ", "; 01118 } 01119 out << vec[vec.dimension()-1]; 01120 } 01121 out << std::string("]"); 01122 return out; 01123 } 01124 template <class T> 01125 std::ostream& operator<<(std::ostream& out, const VectorTranspose<T> &vecT){ 01126 out << vecT.transpose(); 01127 return out; 01128 } 01129 01130 01131 typedef MathVector<int> IntVector; 01132 01133 01134 typedef MathVector<double> RealVector; 01135 01136 01137 // ============================================================================= 01138 // Inversion rules 01139 // ============================================================================= 01140 01141 #define _nemo_vector_inversion_(TYPE) \ 01142 template <> inline MathVector<TYPE> leftAddInverse<MathVector<TYPE> >(const MathVector<TYPE> &value){return -value;} \ 01143 template <> inline MathVector<TYPE> rightAddInverse<MathVector<TYPE> >(const MathVector<TYPE> &value){return -value;} \ 01144 template <> inline MathVector<TYPE> leftSubtractInverse<MathVector<TYPE> >(const MathVector<TYPE> &value){return value;} \ 01145 template <> inline MathVector<TYPE> rightSubtractInverse<MathVector<TYPE> >(const MathVector<TYPE> &value){return -value;} 01146 01147 #ifndef NEMO_PARSED_BY_DOXYGEN 01148 //_nemo_vector_inversion_(float) 01149 //_nemo_vector_inversion_(double) 01150 //_nemo_vector_inversion_(long double) 01151 //_nemo_vector_inversion_(short) 01152 _nemo_vector_inversion_(int) 01153 //_nemo_vector_inversion_(long int) 01154 #endif 01155 01156 #define _nemo_vector_cwise_inversion_(TYPE) \ 01157 _nemo_vector_inversion_(TYPE) \ 01158 template <> inline MathVector<TYPE>::CWise leftMultInverse<>(const MathVector<TYPE>::CWise &value){ \ 01159 return ((MathVector<TYPE>)value).map(((TYPE)1.0)/arg<TYPE>()).cwise(); \ 01160 } \ 01161 template <> inline MathVector<TYPE>::CWise rightMultInverse<>(const MathVector<TYPE>::CWise &value){ \ 01162 return ((MathVector<TYPE>)value).map(((TYPE)1.0)/arg<TYPE>()).cwise(); \ 01163 } \ 01164 template <> inline MathVector<TYPE>::CWise leftDivideInverse<>(const MathVector<TYPE>::CWise &value){ \ 01165 return value; \ 01166 } \ 01167 template <> inline MathVector<TYPE>::CWise rightDivideInverse<>(const MathVector<TYPE>::CWise &value){ \ 01168 return ((MathVector<TYPE>)value).map(((TYPE)1.0)/arg<TYPE>()).cwise(); \ 01169 } 01170 01171 //_nemo_vector_cwise_inversion_(float) 01172 _nemo_vector_cwise_inversion_(double) 01173 //_nemo_vector_cwise_inversion_(long double) 01174 01175 01176 } // end namespace 01177 01178 01179 01180 01181 01182 01183 #endif 01184
Generated on Mon Feb 25 12:49:59 2013 for NemoMath by 1.6.3