OpenWalnut  1.2.5
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
WTensorBase.h
1 //---------------------------------------------------------------------------
2 //
3 // Project: OpenWalnut ( http://www.openwalnut.org )
4 //
5 // Copyright 2009 OpenWalnut Community, BSV@Uni-Leipzig and CNCF@MPI-CBS
6 // For more information see http://www.openwalnut.org/copying
7 //
8 // This file is part of OpenWalnut.
9 //
10 // OpenWalnut is free software: you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as published by
12 // the Free Software Foundation, either version 3 of the License, or
13 // (at your option) any later version.
14 //
15 // OpenWalnut is distributed in the hope that it will be useful,
16 // but WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 // GNU Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public License
21 // along with OpenWalnut. If not, see <http://www.gnu.org/licenses/>.
22 //
23 //---------------------------------------------------------------------------
24 
25 #ifndef WTENSORBASE_H
26 #define WTENSORBASE_H
27 
28 #include <algorithm>
29 #include <map>
30 #include <vector>
31 
32 #include <boost/static_assert.hpp>
33 
34 #include "../WAssert.h"
35 #include "WCompileTimeFunctions.h"
36 #include "WMatrix.h"
37 #include "WValue.h"
38 
39 // TODO(reichenbach): Remove vectors (because of the enum dataSize).
40 
41 // forward declaration of the test classes
42 class WTensorFuncTest;
43 class WTensorBaseTest;
44 class WTensorBaseSymTest;
45 
46 // forward declaration
47 template< std::size_t order, std::size_t dim, typename Data_T >
48 class WTensorBaseSym;
49 
50 // ############################ utility functions ####################################
51 
52 /**
53  * Iterate a position in a multidimensional grid.
54  *
55  * This essentially creates a linear order on all positions in a tensor, where a position
56  * is any valid allocation of indices.
57  *
58  * Example: for a tensor of order 2 and dimension 3, all possible positions are:
59  *
60  * ( 0, 0 ), ( 0, 1 ), ( 0, 2 ), ( 1, 0 ), ( 1, 1 ), ( 1, 2 ), ( 2, 0 ), ( 2, 1 ), ( 2, 2 )
61  *
62  * \param pos The position to be iterated.
63  */
64 template< std::size_t order, std::size_t dim >
65 inline void positionIterateOneStep( std::vector< std::size_t >& pos ) // NOLINT, need a reference here
66 {
67  WAssert( pos.size() >= order, "" );
68 
69  for( std::size_t k = order - 1; k > 0; --k )
70  {
71  if( pos[ k ] == dim - 1)
72  {
73  pos[ k ] = 0;
74  }
75  else
76  {
77  ++( pos[ k ] );
78  return;
79  }
80  }
81  ++( pos[ 0 ] );
82 }
83 
84 /**
85  * Iterate a sorted position in a multidimensional grid.
86  *
87  * This essentially creates a linear order on all sorted positions in a tensor, where a sorted position
88  * is any valid allocation of indices where those indices are in ascending order.
89  *
90  * Example: for a tensor of order 2 and dimension 3, all possible sorted positions are:
91  *
92  * ( 0, 0 ), ( 0, 1 ), ( 0, 2 ), ( 1, 1 ), ( 1, 2 ), ( 2, 2 )
93  *
94  * The number of sorted positions of a tensor matches the number of elements that need to be stored
95  * by the symmetric tensor of the same order and dimension.
96  *
97  * \param pos The sorted(!) position to be iterated.
98  */
99 template< std::size_t order, std::size_t dim >
100 inline void positionIterateSortedOneStep( std::vector< std::size_t >& pos ) // NOLINT, need a reference here
101 {
102  WAssert( pos.size() >= order, "" );
103 
104  for( int k = order - 1; k > -1; --k )
105  {
106  if( pos[ k ] != dim - 1 )
107  {
108  ++( pos[ k ] );
109  for( std::size_t i = k + 1; i < order; ++i )
110  {
111  pos[ i ] = pos[ k ];
112  }
113  return;
114  }
115  }
116 }
117 
118 /**
119  * Same as the version above, using no template arguments.
120  *
121  * \param order The order of the tensor.
122  * \param dim The dimension of the tensor.
123  * \param pos The sorted(!) position to be iterated.
124  */
125 inline void positionIterateSortedOneStep( std::size_t order, std::size_t dim, std::vector< std::size_t >& pos ) // NOLINT, need a reference here
126 {
127  WAssert( pos.size() >= order, "" );
128 
129  for( int k = order - 1; k > -1; --k )
130  {
131  if( pos[ k ] != dim - 1 )
132  {
133  ++( pos[ k ] );
134  for( std::size_t i = k + 1; i < order; ++i )
135  {
136  pos[ i ] = pos[ k ];
137  }
138  return;
139  }
140  }
141 }
142 
143 // ############################# class WTensorBase<> #################################
144 /**
145  * Normal tensor base class.
146  *
147  * \tparam order The order of the tensor.
148  * \tparam dim The dimension of the tensor, i.e. the number of components
149  * in each direction.
150  * \tparam Data_T The datatype of the components, double by default.
151  *
152  * \note The type Data_T may not throw exceptions on construction, destruction or
153  * during any assignment operator.
154  * \note The dimension may not be 0.
155  *
156  * \see WTensor
157  */
158 template< std::size_t order, std::size_t dim, typename Data_T >
160 {
161  /**
162  * For dim == 0, create an artificial compiler error.
163  */
164  BOOST_STATIC_ASSERT( dim != 0 );
165 
166  // Make the appropriate tensor of order + 1 a friend, so it has access to the getPos() member.
167  friend class WTensorBase< order + 1, dim, Data_T >;
168 
169  // make the test class a friend
170  friend class ::WTensorBaseTest;
171 
172  // make the func test class a friend
173  friend class ::WTensorFuncTest;
174 
175 public:
176  /**
177  * Standard constructor.
178  *
179  * All elements are set to Data_T().
180  */
181  WTensorBase();
182 
183  /**
184  * Copy constructor.
185  *
186  * \param t The tensor to copy from.
187  */
188  WTensorBase( WTensorBase const& t ); // NOLINT
189 
190  /**
191  * Copy construct a WTensorBase from a WTensorBaseSym.
192  *
193  * \param t The symmetric tensor to copy from.
194  */
196 
197  /**
198  * Copy operator.
199  *
200  * \param t The tensor to copy from.
201  *
202  * \return *this.
203  */
204  WTensorBase const& operator = ( WTensorBase const& t );
205 
206  /**
207  * Copy operator.
208  *
209  * \param t The symmetric tensor to copy from.
210  *
211  * \return *this.
212  */
214 
215  /**
216  * Get the dimension of this tensor.
217  *
218  * \return The dimension of this tensor.
219  */
220  std::size_t getDimension() const;
221 
222  /**
223  * Get the order of this tensor.
224  *
225  * \return The order of this tensor.
226  */
227  std::size_t getOrder() const;
228 
229  /**
230  * Get the element at a specific position.
231  *
232  * \param indices A std::vector of indices that has a size of at least order.
233  *
234  * \return A reference to the element.
235  */
236  template< typename Index_T >
237  Data_T& operator[] ( std::vector< Index_T > const& indices );
238 
239  /**
240  * Get the element at a specific position.
241  *
242  * \param indices A std::vector of indices that has a size of at least order.
243  *
244  * \return A reference to the element.
245  */
246  template< typename Index_T >
247  Data_T const& operator[] ( std::vector< Index_T > const& indices ) const;
248 
249  /**
250  * Get the element at a specific position.
251  *
252  * \param indices An array of indices that has a size of at least order.
253  *
254  * \return A reference to the element.
255  *
256  * \note The array must have a length of at least order.
257  */
258  template< typename Index_T >
259  Data_T& operator[] ( Index_T indices[] );
260 
261  /**
262  * Get the element at a specific position.
263  *
264  * \param indices An array of indices that has a size of at least order.
265  *
266  * \return A reference to the element.
267  *
268  * \note The array must have a length of at least order.
269  */
270  template< typename Index_T >
271  Data_T const& operator[] ( Index_T indices[] ) const;
272 
273  /**
274  * Compare this WTensorBase to another one.
275  *
276  * \param other The WBensorBase to compare to.
277  *
278  * \return True, iff this tensors' elements are equal to another tensors' elements.
279  */
280  bool operator == ( WTensorBase const& other ) const;
281 
282  /**
283  * Compare this WTensorBase to another one.
284  *
285  * \param other The WBensorBase to compare to.
286  *
287  * \return True, iff this tensors' elements are not equal to another tensors' elements.
288  */
289  bool operator != ( WTensorBase const& other ) const;
290 
291 private:
292  /**
293  * Calculate the position of the element in the data vector. The function
294  * is build recursively at compile-time.
295  *
296  * \param pos An array of indices.
297  *
298  * \return The position of the element.
299  */
300  template< typename Index_T >
301  static inline std::size_t getPos( Index_T pos[] );
302 
303  /**
304  * Stores all elements.
305  */
306  std::vector< Data_T > m_data;
307 
308  /**
309  * Declare a compile-time constant as enum and not as static constant.
310  */
311  enum
312  {
313  /**
314  * The number of elements to store.
315  */
317  };
318 };
319 
320 template< std::size_t order, std::size_t dim, typename Data_T >
322  : m_data( dataSize, Data_T() )
323 {
324 }
325 
326 template< std::size_t order, std::size_t dim, typename Data_T >
328  : m_data( t.m_data )
329 {
330 }
331 
332 template< std::size_t order, std::size_t dim, typename Data_T >
334  : m_data( dataSize )
335 {
336  *this = t;
337 }
338 
339 template< std::size_t order, std::size_t dim, typename Data_T >
341 {
342  m_data = t.m_data;
343  return *this;
344 }
345 
346 template< std::size_t order, std::size_t dim, typename Data_T >
348 {
349  std::vector< std::size_t > pos( order, 0 );
350 
351  for( std::size_t k = 0; k < dataSize; ++k )
352  {
353  ( *this )[ pos ] = t[ pos ];
354  positionIterateOneStep< order, dim >( pos );
355  }
356 
357  return *this;
358 }
359 
360 template< std::size_t order, std::size_t dim, typename Data_T >
362 {
363  return dim;
364 }
365 
366 template< std::size_t order, std::size_t dim, typename Data_T >
368 {
369  return order;
370 }
371 
372 template< std::size_t order, std::size_t dim, typename Data_T >
373 template< typename Index_T >
374 std::size_t WTensorBase< order, dim, Data_T >::getPos( Index_T pos[] )
375 {
376  return WTensorBase< order - 1, dim, Data_T >::getPos( pos ) * dim + static_cast< std::size_t >( pos[ order - 1 ] );
377 }
378 
379 template< std::size_t order, std::size_t dim, typename Data_T >
380 template< typename Index_T >
381 Data_T& WTensorBase< order, dim, Data_T >::operator[] ( std::vector< Index_T > const& indices )
382 {
383  return const_cast< Data_T& >( static_cast< WTensorBase const& >( *this ).operator[] ( indices ) );
384 }
385 
386 template< std::size_t order, std::size_t dim, typename Data_T >
387 template< typename Index_T >
388 Data_T const& WTensorBase< order, dim, Data_T >::operator[] ( std::vector< Index_T > const& indices ) const
389 {
390  WAssert( indices.size() >= order, "" );
391  return operator[] ( &indices[ 0 ] );
392 }
393 
394 template< std::size_t order, std::size_t dim, typename Data_T >
395 template< typename Index_T >
397 {
398  return const_cast< Data_T& >( static_cast< WTensorBase const& >( *this ).operator[] ( indices ) );
399 }
400 
401 template< std::size_t order, std::size_t dim, typename Data_T >
402 template< typename Index_T >
403 Data_T const& WTensorBase< order, dim, Data_T >::operator[] ( Index_T indices[] ) const
404 {
405  for( std::size_t k = 0; k < order; ++k )
406  {
407  WAssert( static_cast< std::size_t >( indices[ k ] ) < dim, "" );
408  }
409  std::size_t p = getPos( indices );
410  return m_data[ p ];
411 }
412 
413 template< std::size_t order, std::size_t dim, typename Data_T >
415 {
416  return m_data == other.m_data;
417 }
418 
419 template< std::size_t order, std::size_t dim, typename Data_T >
421 {
422  return m_data != other.m_data;
423 }
424 
425 // ######################### WTensorBase for order == 0 ##########################
426 
427 /**
428  * \tparam dim The dimension of the tensor, i.e. the number of components
429  * in each direction.
430  * \tparam Data_T The datatype of the components, double by default.
431 
432  * Specialization for order = 0. This essentially encapsulates a scalar. The purpose of
433  * this specialization is compatibility for generic tensor functions.
434  */
435 template< std::size_t dim, typename Data_T >
436 class WTensorBase< 0, dim, Data_T >
437 {
438  // Make the appropriate tensor of order + 1 a friend, so it has access to the getPos() member.
439  friend class WTensorBase< 1, dim, Data_T >;
440 
441  // make the test class a friend
442  friend class ::WTensorBaseTest;
443 
444  // make the func test class a friend
445  friend class ::WTensorFuncTest;
446 
447 public:
448  // implementations are inline as they are trivial
449  /**
450  * Standard constructor.
451  *
452  * All elements are set to Data_T().
453  */
455  : m_data( Data_T() )
456  {
457  }
458 
459  /**
460  * Copy constructor.
461  *
462  * \param t The tensor to copy from.
463  */
464  WTensorBase( WTensorBase const& t ) // NOLINT
465  : m_data( t.m_data )
466  {
467  }
468 
469  /**
470  * Copy construct a WTensorBase from a WTensorBaseSym.
471  *
472  * \param t The symmetric tensor to copy from.
473  */
475  : m_data()
476  {
477  m_data = t.operator[]< std::size_t >( NULL );
478  }
479 
480  /**
481  * Copy operator.
482  *
483  * \param t The tensor to copy from.
484  *
485  * \return *this.
486  */
488  {
489  m_data = t.m_data;
490  return *this;
491  }
492 
493  /**
494  * Copy operator.
495  *
496  * \param t The symmetric tensor to copy from.
497  *
498  * \return *this.
499  */
501  {
502  m_data = t.operator[]< std::size_t >( NULL );
503  return *this;
504  }
505 
506  /**
507  * Get the dimension of this tensor.
508  *
509  * \return The dimension of this tensor.
510  */
511  std::size_t getDimension() const
512  {
513  return dim;
514  }
515 
516  /**
517  * Get the order of this tensor.
518  *
519  * \return The order of this tensor.
520  */
521  std::size_t getOrder() const
522  {
523  return 0;
524  }
525 
526  /**
527  * Get the value of this scalar.
528  *
529  * \return A reference to the element.
530  */
531  template< typename Index_T >
532  Data_T& operator[] ( std::vector< Index_T > const& /* indices */ )
533  {
534  return m_data;
535  }
536 
537  /**
538  * Get the value of this scalar.
539  *
540  * \return A reference to the element.
541  */
542  template< typename Index_T >
543  Data_T const& operator[] ( std::vector< Index_T > const& /* indices */ ) const
544  {
545  return m_data;
546  }
547 
548  /**
549  * Get the value of this scalar.
550  *
551  * \return A reference to the element.
552  */
553  template< typename Index_T >
554  Data_T& operator[] ( Index_T[] /* indices */ )
555  {
556  return m_data;
557  }
558 
559  /**
560  * Get the value of this scalar.
561  *
562  * \return A reference to the element.
563  */
564  template< typename Index_T >
565  Data_T const& operator[] ( Index_T[] /* indices */ ) const
566  {
567  return m_data;
568  }
569 
570  /**
571  * Compare this WTensorBase to another one.
572  *
573  * \param other The WBensorBase to compare to.
574  *
575  * \return True, iff this tensors' elements are equal to another tensors' elements.
576  */
577  bool operator == ( WTensorBase const& other ) const
578  {
579  return m_data == other.m_data;
580  }
581 
582  /**
583  * Compare this WTensorBase to another one.
584  *
585  * \param other The WBensorBase to compare to.
586  *
587  * \return True, iff this tensors' elements are not equal to another tensors' elements.
588  */
589  bool operator != ( WTensorBase const& other ) const
590  {
591  return m_data != other.m_data;
592  }
593 
594 private:
595  /**
596  * Calculate the position of the element in the data vector. This
597  * is essentially the standard case of the recursion.
598  *
599  * \return 0.
600  */
601  template< typename Index_T >
602  static inline std::size_t getPos( Index_T[] /* pos */ )
603  {
604  return 0;
605  }
606 
607  /**
608  * Stores the value.
609  */
610  Data_T m_data;
611 
612  /**
613  * Declare a compile-time constant as enum and not as static constant.
614  */
615  enum
616  {
617  /**
618  * The number of elements to store.
619  */
621  };
622 };
623 
624 // ################################# class WTensorBaseSym<> #####################################
625 
626 /**
627  * Symmetric tensor base class.
628  *
629  * \tparam order The order of the tensor.
630  * \tparam dim The dimension of the tensor, i.e. the number of components
631  * in each direction.
632  * \tparam Data_T The datatype of the components, double by default.
633  *
634  * \note The type Data_T may not throw exceptions on construction, destruction or
635  * during any assignment operator.
636  * \note The dimension may not be 0.
637  *
638  * \see WTensorSym
639  */
640 template< std::size_t order, std::size_t dim, typename Data_T >
642 {
643  /**
644  * For dim == 0, create an artificial compiler error.
645  */
646  BOOST_STATIC_ASSERT( dim != 0 );
647 
648  // make the test class a friend
649  friend class ::WTensorBaseSymTest;
650 
651  // make the func test class a friend
652  friend class ::WTensorFuncTest;
653 
654 public:
655  /**
656  * Standard constructor.
657  *
658  * All elements are set to Data_T().
659  */
660  WTensorBaseSym();
661 
662  /**
663  * Constructs the symmetrical tensor and initialize with the given data.
664  *
665  * \param data The components of the symmetrical tensor: Take care of the
666  * ordering of the components to match the ordering in \see m_data.
667  */
668  explicit WTensorBaseSym( const WValue< Data_T >& data );
669 
670  /**
671  * Copy constructor.
672  *
673  * \param t The tensor to copy from.
674  */
675  WTensorBaseSym( WTensorBaseSym const& t ); // NOLINT
676 
677  /**
678  * Copy operator.
679  *
680  * \param t The tensor to copy from.
681  *
682  * \return *this.
683  */
684  WTensorBaseSym const& operator = ( WTensorBaseSym const& t );
685 
686  /**
687  * Get the dimension of this tensor.
688  *
689  * \return The dimension of this tensor.
690  */
691  std::size_t getDimension() const;
692 
693  /**
694  * Get the order of this tensor.
695  *
696  * \return The order of this tensor.
697  */
698  std::size_t getOrder() const;
699 
700  /**
701  * Get the element at a specific position.
702  *
703  * \param indices A std::vector of indices that has a size of at least order.
704  *
705  * \return A reference to the element.
706  */
707  template< typename Index_T >
708  Data_T& operator[] ( std::vector< Index_T > const& indices );
709 
710  /**
711  * Get the element at a specific position.
712  *
713  * \param indices A std::vector of indices that has a size of at least order.
714  *
715  * \return A reference to the element.
716  */
717  template< typename Index_T >
718  Data_T const& operator[] ( std::vector< Index_T > const& indices ) const;
719 
720  /**
721  * Get the element at a specific position.
722  *
723  * \param indices An array of indices that has a size of at least order.
724  *
725  * \return A reference to the element.
726  *
727  * \note No bounds checking is performed.
728  */
729  template< typename Index_T >
730  Data_T& operator[] ( Index_T indices[] );
731 
732  /**
733  * Get the element at a specific position.
734  *
735  * \param indices An array of indices that has a size of at least order.
736  *
737  * \return A reference to the element.
738  *
739  * \note No bounds checking is performed.
740  */
741  template< typename Index_T >
742  Data_T const& operator[] ( Index_T indices[] ) const;
743 
744  /**
745  * Compare this WTensorBaseSym to another one.
746  *
747  * \param other The WTensorBaseSym to compare to.
748  *
749  * \return True, iff this tensors' elements are equal to another tensors' elements.
750  */
751  bool operator == ( WTensorBaseSym const& other ) const;
752 
753  /**
754  * Compare this WTensorBaseSym to another one.
755  *
756  * \param other The WTensorBaseSym to compare to.
757  *
758  * \return True, iff this tensors' elements are not equal to another tensors' elements.
759  */
760  bool operator != ( WTensorBaseSym const& other ) const;
761 
762 private:
763  /**
764  * Stores the elements of this tensor lexicographical ordered on their
765  * indices, where for each set of permutations the lexicographical lowest
766  * index is used.
767  */
768  std::vector< Data_T > m_data;
769 
770  /**
771  * Declare a compile-time constant as enum and not as static constant.
772  */
773  enum
774  {
775  /**
776  * The number of elements to store.
777  */
778  dataSize = WBinom< order + dim - 1, order >::value
779  };
780 
781  /**
782  * A class that maps symmetric tensor indices to vector positions.
783  */
785  {
786  public:
787  /**
788  * Standard constructor. The mapping is calculated here.
789  */
790  PositionIndexer();
791 
792  /**
793  * Get the mapped position.
794  *
795  * \param pos An array of indices.
796  *
797  * \return The position that corresponds to the indices.
798  */
799  template< typename Index_T >
800  inline std::size_t operator[] ( Index_T pos[] ) const;
801 
802  /**
803  * Declare a compile-time constant as enum and not as static constant.
804  */
805  enum
806  {
807  /**
808  * The number of data elements.
809  */
810  dataSize = WBinom< order + dim - 1, order >::value
811  };
812 
813  private:
814  /**
815  * Maps the indices to a vector element position.
816  */
818  };
819 
820  /**
821  * A static PositionIndexer that maps tensor indices to vector positions.
822  */
823  static PositionIndexer const m_indexer;
824 };
825 
826 // initialize the indexer object as a static object
827 template< std::size_t order, std::size_t dim, typename Data_T >
829 
830 // ######################## impl of WTensorBaseSym::PositionIndexer #####################################
831 
832 template< std::size_t order, std::size_t dim, typename Data_T >
834 {
835  // the map uses lexical ordering of vectors
836  std::map< std::vector< std::size_t >, std::size_t > m;
837 
838  // fill the map with all possible combinations of indices, where
839  // every combination of indices appears in ascending order of indices
840  std::vector< std::size_t > pos( order, 0 );
841  for( std::size_t k = 0; k < dataSize; ++k )
842  {
843  // enumerate the position
844  m[ pos ] = k;
845 
846  // get the next sorted combination
847  positionIterateSortedOneStep< order, dim >( pos );
848  }
849 
850  // now iterate all possible sets of indices
851  pos = std::vector< std::size_t >( order, 0 );
852  std::vector< std::size_t > _p( order, 0 );
853  for( std::size_t k = 0; k < WPower< dim, order >::value; ++k )
854  {
855  _p = pos;
856 
857  // sort the indices in _p
858  std::sort( _p.begin(), _p.end() );
859 
860  // now map the arbitrary ordered indices to the position of the ordered set in m (and thus in m_data)
861  m_positions[ pos ] = m[ _p ];
862 
863  // the map should already know the sorted position,
864  // it should never be added by std::map::operator [] at this point
865  WAssert( m.size() == dataSize, "" );
866 
867  // get the next position
868  positionIterateOneStep< order, dim >( pos );
869  }
870 }
871 
872 template< std::size_t order, std::size_t dim, typename Data_T >
873 template< typename Index_T >
875 {
876  return m_positions[ pos ];
877 }
878 
879 // ######################## impl of WTensorBaseSym #####################################
880 
881 template< std::size_t order, std::size_t dim, typename Data_T >
883  : m_data( dataSize, Data_T() )
884 {
885 }
886 
887 template< std::size_t order, std::size_t dim, typename Data_T >
889  : m_data( &data[0], &data[0] + data.size() )
890 {
891  WAssert( dataSize == m_data.size(), "Number of given components does not match the order and dimension of this symmetric tensor" );
892 }
893 
894 template< std::size_t order, std::size_t dim, typename Data_T >
896  : m_data( t.m_data )
897 {
898 }
899 
900 template< std::size_t order, std::size_t dim, typename Data_T >
902 {
903  m_data = t.m_data;
904  return *this;
905 }
906 
907 template< std::size_t order, std::size_t dim, typename Data_T >
909 {
910  return dim;
911 }
912 
913 template< std::size_t order, std::size_t dim, typename Data_T >
915 {
916  return order;
917 }
918 
919 template< std::size_t order, std::size_t dim, typename Data_T >
920 template< typename Index_T >
921 Data_T& WTensorBaseSym< order, dim, Data_T >::operator[] ( std::vector< Index_T > const& indices )
922 {
923  return const_cast< Data_T& >( static_cast< WTensorBaseSym const& >( *this ).operator[] ( indices ) );
924 }
925 
926 template< std::size_t order, std::size_t dim, typename Data_T >
927 template< typename Index_T >
928 Data_T const& WTensorBaseSym< order, dim, Data_T >::operator[] ( std::vector< Index_T > const& indices ) const
929 {
930  WAssert( indices.size() >= order, "" );
931  return operator[] ( &indices[ 0 ] );
932 }
933 
934 template< std::size_t order, std::size_t dim, typename Data_T >
935 template< typename Index_T >
937 {
938  return const_cast< Data_T& >( static_cast< WTensorBaseSym const& >( *this ).operator[] ( indices ) );
939 }
940 
941 template< std::size_t order, std::size_t dim, typename Data_T >
942 template< typename Index_T >
943 Data_T const& WTensorBaseSym< order, dim, Data_T >::operator[] ( Index_T indices[] ) const
944 {
945  std::size_t p = m_indexer[ indices ];
946  WAssert( p < m_data.size(), "" );
947  return m_data[ p ];
948 }
949 
950 template< std::size_t order, std::size_t dim, typename Data_T >
952 {
953  return m_data == other.m_data;
954 }
955 
956 template< std::size_t order, std::size_t dim, typename Data_T >
958 {
959  return m_data != other.m_data;
960 }
961 
962 // ######################### WTensorBaseSym for order == 0 ##########################
963 
964 /**
965  * \tparam dim The dimension of the tensor, i.e. the number of components
966  * in each direction.
967  * \tparam Data_T The datatype of the components, double by default.
968 
969  * Specialization for order = 0. This essentially encapsulates a scalar. The purpose of
970  * this specialization is compatibility for generic tensor functions.
971  */
972 template< std::size_t dim, typename Data_T >
973 class WTensorBaseSym< 0, dim, Data_T >
974 {
975  friend class ::WTensorBaseSymTest;
976  friend class ::WTensorFuncTest;
977 
978 public:
979  // implementations are inline as they are trivial
980  /**
981  * Standard constructor.
982  *
983  * All elements are set to Data_T().
984  */
986  : m_data( Data_T() )
987  {
988  }
989 
990  /**
991  * Copy constructor.
992  *
993  * \param t The tensor to copy from.
994  */
995  WTensorBaseSym( WTensorBaseSym const& t ) // NOLINT
996  : m_data( t.m_data )
997  {
998  }
999 
1000  /**
1001  * Copy operator.
1002  *
1003  * \param t The tensor to copy from.
1004  *
1005  * \return *this.
1006  */
1008  {
1009  m_data = t.m_data;
1010  return *this;
1011  }
1012 
1013  /**
1014  * Get the dimension of this tensor.
1015  *
1016  * \return The dimension of this tensor.
1017  */
1018  std::size_t getDimension() const
1019  {
1020  return dim;
1021  }
1022 
1023  /**
1024  * Get the order of this tensor.
1025  *
1026  * \return The order of this tensor.
1027  */
1028  std::size_t getOrder() const
1029  {
1030  return 0;
1031  }
1032 
1033  /**
1034  * Get the value of this scalar.
1035  *
1036  * \return A reference to the element.
1037  */
1038  template< typename Index_T >
1039  Data_T& operator[] ( std::vector< Index_T > const& /* indices */ )
1040  {
1041  return m_data;
1042  }
1043 
1044  /**
1045  * Get the value of this scalar.
1046  *
1047  * \return A reference to the element.
1048  */
1049  template< typename Index_T >
1050  Data_T const& operator[] ( std::vector< Index_T > const& /* indices */ ) const
1051  {
1052  return m_data;
1053  }
1054 
1055  /**
1056  * Get the value of this scalar.
1057  *
1058  * \return A reference to the element.
1059  */
1060  template< typename Index_T >
1061  Data_T& operator[] ( Index_T[] /* indices */ )
1062  {
1063  return m_data;
1064  }
1065 
1066  /**
1067  * Get the value of this scalar.
1068  *
1069  * \return A reference to the element.
1070  */
1071  template< typename Index_T >
1072  Data_T const& operator[] ( Index_T[] /* indices */ ) const
1073  {
1074  return m_data;
1075  }
1076 
1077  /**
1078  * Compare this WTensorBase to another one.
1079  *
1080  * \param other The WBensorBase to compare to.
1081  *
1082  * \return True, iff this tensors' elements are equal to another tensors' elements.
1083  */
1084  bool operator == ( WTensorBaseSym const& other ) const
1085  {
1086  return m_data == other.m_data;
1087  }
1088 
1089  /**
1090  * Compare this WTensorBase to another one.
1091  *
1092  * \param other The WBensorBase to compare to.
1093  *
1094  * \return True, iff this tensors' elements are not equal to another tensors' elements.
1095  */
1096  bool operator != ( WTensorBaseSym const& other ) const
1097  {
1098  return m_data != other.m_data;
1099  }
1100 
1101 private:
1102  /**
1103  * Stores the value.
1104  */
1105  Data_T m_data;
1106 
1107  /**
1108  * Declare a compile-time constant as enum and not as static constant.
1109  */
1110  enum
1111  {
1112  /**
1113  * The number of elements to store.
1114  */
1116  };
1117 };
1118 
1119 // ################################### class WTensorFunc<> ######################################
1120 
1121 /**
1122  * Implements functions that should only be defined for certain values of order.
1123  *
1124  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1125  * \tparam order The order of the tensor.
1126  * \tparam dim The dimension of the tensor, i.e. the number of components
1127  * in each direction.
1128  * \tparam Data_T The datatype of the components, double by default.
1129  */
1130 //next line is nolint because brainlint cannot find the declaration of TensorBase_T
1131 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t order, std::size_t dim, typename Data_T > //NOLINT
1132 class WTensorFunc : public TensorBase_T< order, dim, Data_T >
1133 {
1134 };
1135 
1136 /**
1137  * Implements the operator () for an order of 6.
1138  *
1139  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1140  * \tparam dim The dimension of the tensor, i.e. the number of components
1141  * in each direction.
1142  * \tparam Data_T The datatype of the components, double by default.
1143  */
1144 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1145 class WTensorFunc< TensorBase_T, 6, dim, Data_T > : public TensorBase_T< 6, dim, Data_T >
1146 {
1147 public:
1148  /**
1149  * Access operator.
1150  *
1151  * \param i0 An index.
1152  * \param i1 An index.
1153  * \param i2 An index.
1154  * \param i3 An index.
1155  * \param i4 An index.
1156  * \param i5 An index.
1157  *
1158  * \return A reference to the element.
1159  */
1160  Data_T& operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3, std::size_t i4, std::size_t i5 );
1161 
1162  /**
1163  * Access operator.
1164  *
1165  * \param i0 An index.
1166  * \param i1 An index.
1167  * \param i2 An index.
1168  * \param i3 An index.
1169  * \param i4 An index.
1170  * \param i5 An index.
1171  *
1172  * \return A reference to the element.
1173  */
1174  Data_T const& operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3, std::size_t i4, std::size_t i5 ) const;
1175 };
1176 
1177 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1178 Data_T& WTensorFunc< TensorBase_T, 6, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2,
1179  std::size_t i3, std::size_t i4, std::size_t i5 )
1180 {
1181  return const_cast< Data_T& >( static_cast< WTensorFunc const& >( *this ).operator() ( i0, i1, i2, i3, i4, i5 ) );
1182 }
1183 
1184 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1185 Data_T const& WTensorFunc< TensorBase_T, 6, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2,
1186  std::size_t i3, std::size_t i4, std::size_t i5 ) const
1187 {
1188  std::size_t p[] = { i0, i1, i2, i3, i4, i5 };
1190 }
1191 
1192 
1193 
1194 
1195 /**
1196  * Implements the operator () for an order of 5.
1197  *
1198  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1199  * \tparam dim The dimension of the tensor, i.e. the number of components
1200  * in each direction.
1201  * \tparam Data_T The datatype of the components, double by default.
1202  */
1203 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1204 class WTensorFunc< TensorBase_T, 5, dim, Data_T > : public TensorBase_T< 5, dim, Data_T >
1205 {
1206 public:
1207  /**
1208  * Access operator.
1209  *
1210  * \param i0 An index.
1211  * \param i1 An index.
1212  * \param i2 An index.
1213  * \param i3 An index.
1214  * \param i4 An index.
1215  *
1216  * \return A reference to the element.
1217  */
1218  Data_T& operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3, std::size_t i4 );
1219 
1220  /**
1221  * Access operator.
1222  *
1223  * \param i0 An index.
1224  * \param i1 An index.
1225  * \param i2 An index.
1226  * \param i3 An index.
1227  * \param i4 An index.
1228  *
1229  * \return A reference to the element.
1230  */
1231  Data_T const& operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3, std::size_t i4 ) const;
1232 };
1233 
1234 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1235 Data_T& WTensorFunc< TensorBase_T, 5, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2,
1236  std::size_t i3, std::size_t i4 )
1237 {
1238  return const_cast< Data_T& >( static_cast< WTensorFunc const& >( *this ).operator() ( i0, i1, i2, i3, i4 ) );
1239 }
1240 
1241 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1242 Data_T const& WTensorFunc< TensorBase_T, 5, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2,
1243  std::size_t i3, std::size_t i4 ) const
1244 {
1245  std::size_t p[] = { i0, i1, i2, i3, i4 };
1247 }
1248 
1249 
1250 
1251 
1252 /**
1253  * Implements the operator () for an order of 4.
1254  *
1255  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1256  * \tparam dim The dimension of the tensor, i.e. the number of components
1257  * in each direction.
1258  * \tparam Data_T The datatype of the components, double by default.
1259  */
1260 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1261 class WTensorFunc< TensorBase_T, 4, dim, Data_T > : public TensorBase_T< 4, dim, Data_T >
1262 {
1263 public:
1264  /**
1265  * Access operator.
1266  *
1267  * \param i0 An index.
1268  * \param i1 An index.
1269  * \param i2 An index.
1270  * \param i3 An index.
1271  *
1272  * \return A reference to the element.
1273  */
1274  Data_T& operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3 );
1275 
1276  /**
1277  * Access operator.
1278  *
1279  * \param i0 An index.
1280  * \param i1 An index.
1281  * \param i2 An index.
1282  * \param i3 An index.
1283  *
1284  * \return A reference to the element.
1285  */
1286  Data_T const& operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3 ) const;
1287 };
1288 
1289 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1290 Data_T& WTensorFunc< TensorBase_T, 4, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3 )
1291 {
1292  return const_cast< Data_T& >( static_cast< WTensorFunc const& >( *this ).operator() ( i0, i1, i2, i3 ) );
1293 }
1294 
1295 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1296 Data_T const& WTensorFunc< TensorBase_T, 4, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2, std::size_t i3 ) const
1297 {
1298  std::size_t p[] = { i0, i1, i2, i3 };
1300 }
1301 
1302 
1303 
1304 /**
1305  * Implements the operator () for an order of 3.
1306  *
1307  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1308  * \tparam dim The dimension of the tensor, i.e. the number of components
1309  * in each direction.
1310  * \tparam Data_T The datatype of the components, double by default.
1311  */
1312 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1313 class WTensorFunc< TensorBase_T, 3, dim, Data_T > : public TensorBase_T< 3, dim, Data_T >
1314 {
1315 public:
1316  /**
1317  * Access operator.
1318  *
1319  * \param i0 An index.
1320  * \param i1 An index.
1321  * \param i2 An index.
1322  *
1323  * \return A reference to the element.
1324  */
1325  Data_T& operator() ( std::size_t i0, std::size_t i1, std::size_t i2 );
1326 
1327  /**
1328  * Access operator.
1329  *
1330  * \param i0 An index.
1331  * \param i1 An index.
1332  * \param i2 An index.
1333  *
1334  * \return A reference to the element.
1335  */
1336  Data_T const& operator() ( std::size_t i0, std::size_t i1, std::size_t i2 ) const;
1337 };
1338 
1339 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1340 Data_T& WTensorFunc< TensorBase_T, 3, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2 )
1341 {
1342  return const_cast< Data_T& >( static_cast< WTensorFunc const& >( *this ).operator() ( i0, i1, i2 ) );
1343 }
1344 
1345 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1346 Data_T const& WTensorFunc< TensorBase_T, 3, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1, std::size_t i2 ) const
1347 {
1348  std::size_t p[] = { i0, i1, i2 };
1350 }
1351 
1352 /**
1353  * Implements the operator () for an order of 2 as well as casts to WMatrix.
1354  *
1355  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1356  * \tparam dim The dimension of the tensor, i.e. the number of components
1357  * in each direction.
1358  * \tparam Data_T The datatype of the components, double by default.
1359  */
1360 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1361 class WTensorFunc< TensorBase_T, 2, dim, Data_T > : public TensorBase_T< 2, dim, Data_T >
1362 {
1363 public:
1364  /**
1365  * Default constructor.
1366  */
1367  WTensorFunc();
1368 
1369  /**
1370  * Initializes the tensor with the given data.
1371  *
1372  * \param data Components in same ordering as the components of the TensorBase class.
1373  */
1374  explicit WTensorFunc( const WValue< Data_T >& data );
1375 
1376  /**
1377  * Access operator.
1378  *
1379  * \param i0 An index.
1380  * \param i1 An index.
1381  *
1382  * \return A reference to the element.
1383  */
1384  Data_T& operator() ( std::size_t i0, std::size_t i1 );
1385 
1386  /**
1387  * Access operator.
1388  *
1389  * \param i0 An index.
1390  * \param i1 An index.
1391  *
1392  * \return A reference to the element.
1393  */
1394  Data_T const& operator() ( std::size_t i0, std::size_t i1 ) const;
1395 
1396  /**
1397  * Cast this 2nd-order tensor into a WMatrix.
1398  */
1399  operator WMatrix< Data_T >() const;
1400 };
1401 
1402 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T >
1404  : TensorBase_T< 2, dim, Data_T >()
1405 {
1406 }
1407 
1408 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T >
1410  : TensorBase_T< 2, dim, Data_T >( data )
1411 {
1412 }
1413 
1414 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1415 Data_T& WTensorFunc< TensorBase_T, 2, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1 )
1416 {
1417  return const_cast< Data_T& >( static_cast< WTensorFunc const& >( *this ).operator() ( i0, i1 ) );
1418 }
1419 
1420 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1421 Data_T const& WTensorFunc< TensorBase_T, 2, dim, Data_T >::operator() ( std::size_t i0, std::size_t i1 ) const
1422 {
1423  std::size_t p[] = { i0, i1 };
1425 }
1426 
1427 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1429 {
1430  WMatrix< Data_T > m( dim, dim );
1431  for( std::size_t i = 0; i < dim; ++i )
1432  {
1433  for( std::size_t j = 0; j < dim; ++j )
1434  {
1435  m( i, j ) = this->operator() ( i, j );
1436  }
1437  }
1438  return m;
1439 }
1440 
1441 /**
1442  * Implements the operator () for an order of 1 as well as a cast to WValue.
1443  *
1444  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1445  * \tparam dim The dimension of the tensor, i.e. the number of components
1446  * in each direction.
1447  * \tparam Data_T The datatype of the components, double by default.
1448  */
1449 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1450 class WTensorFunc< TensorBase_T, 1, dim, Data_T > : public TensorBase_T< 1, dim, Data_T >
1451 {
1452 public:
1453  /**
1454  * Access operator.
1455  *
1456  * \param i0 An index.
1457  *
1458  * \return A reference to the element.
1459  */
1460  Data_T& operator() ( std::size_t i0 );
1461 
1462  /**
1463  * Access operator.
1464  *
1465  * \param i0 An index.
1466  *
1467  * \return A reference to the element.
1468  */
1469  Data_T const& operator() ( std::size_t i0 ) const;
1470 
1471  /**
1472  * Cast this tensor to a WValue.
1473  */
1474  operator WValue< Data_T > () const;
1475 };
1476 
1477 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1479 {
1480  return const_cast< Data_T& >( static_cast< WTensorFunc const& >( *this ).operator() ( i0 ) );
1481 }
1482 
1483 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1485 {
1486  std::size_t p[] = { i0 };
1488 }
1489 
1490 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1492 {
1493  WValue< Data_T > v( dim );
1494 
1495  for( std::size_t k = 0; k < dim; ++k )
1496  {
1497  v[ k ] = this->operator() ( k );
1498  }
1499 
1500  return v;
1501 }
1502 
1503 /**
1504  * Implements the operator () for an order of 0 and a cast to Data_T.
1505  *
1506  * \tparam TensorBase_T Either WTensorBase<> or WTensorBaseSym<>
1507  * \tparam dim The dimension of the tensor, i.e. the number of components
1508  * in each direction.
1509  * \tparam Data_T The datatype of the components, double by default.
1510  */
1511 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1512 class WTensorFunc< TensorBase_T, 0, dim, Data_T > : public TensorBase_T< 0, dim, Data_T >
1513 {
1514 public:
1515  /**
1516  * Access operator.
1517  *
1518  * \return A reference to the element.
1519  */
1520  Data_T& operator() ();
1521 
1522  /**
1523  * Access operator.
1524  *
1525  * \return A reference to the element.
1526  */
1527  Data_T const& operator() () const;
1528 
1529  /**
1530  * Cast this tensor to a Data_T.
1531  */
1532  operator Data_T() const;
1533 };
1534 
1535 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1537 {
1538  return const_cast< Data_T& >( static_cast< WTensorFunc const& >( *this ).operator() () );
1539 }
1540 
1541 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1543 {
1544 #ifdef _MSC_VER
1545  return TensorBase_T< 0, dim, Data_T >::operator[]< std::size_t >( NULL );
1546 #else
1547  return TensorBase_T< 0, dim, Data_T >::template operator[]< std::size_t >( NULL );
1548 #endif
1549 }
1550 
1551 template< template< std::size_t, std::size_t, typename > class TensorBase_T, std::size_t dim, typename Data_T > //NOLINT
1553 {
1554 #ifdef _MSC_VER
1555  return TensorBase_T< 0, dim, Data_T >::operator[]< std::size_t >( NULL );
1556 #else
1557  return TensorBase_T< 0, dim, Data_T >::template operator[]< std::size_t >( NULL );
1558 #endif
1559 }
1560 
1561 #endif // WTENSORBASE_H