1 /*
2   This file is part of MADNESS.
3 
4   Copyright (C) 2007,2010 Oak Ridge National Laboratory
5 
6   This program is free software; you can redistribute it and/or modify
7   it under the terms of the GNU General Public License as published by
8   the Free Software Foundation; either version 2 of the License, or
9   (at your option) any later version.
10 
11   This program is distributed in the hope that it will be useful,
12   but WITHOUT ANY WARRANTY; without even the implied warranty of
13   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14   GNU General Public License for more details.
15 
16   You should have received a copy of the GNU General Public License
17   along with this program; if not, write to the Free Software
18   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 
20   For more information please contact:
21 
22   Robert J. Harrison
23   Oak Ridge National Laboratory
24   One Bethel Valley Road
25   P.O. Box 2008, MS-6367
26 
27   email: harrisonrj@ornl.gov
28   tel:   865-241-3937
29   fax:   865-572-0680
30 
31 
32   $Id$
33 */
34 
35 
36 #include <iostream>
37 
38 #include <algorithm>
39 
40 #include "basetensor.h"
41 
42 #include "tensorexcept.h"
43 
44 /// \file basetensor.cc
45 /// \brief Implements BaseTensor
46 
47 namespace madness {
48 
49 #ifdef TENSOR_INSTANCE_COUNT
50   madness::AtomicInt BaseTensor::instance_count;
51 #endif
52 
53 
54     /// Reshape the size and number of dimensions.
55 
56     /// Modifies the current tensor to have the number and size of
57     /// dimensions as described in the \c vector \c d .  The total number
58     /// of elements must be the same before and after, and the current
59     /// tensor must be contiguous.
reshape_inplace(int nd,const long * d)60     void BaseTensor::reshape_inplace(int nd, const long* d) {
61         TENSOR_ASSERT(iscontiguous(),
62                       "cannot reshape non-contiguous tensor ... consider fuse/splitdim",
63                       0,this);
64         long newsize=1;
65         for (long i=0; i<nd; ++i) newsize *= d[i];
66         TENSOR_ASSERT(_size == newsize,"old and new sizes do not match",_size,this);
67         set_dims_and_size(nd,&(d[0]));
68     }
69 
70     /// Reshape the size and number of dimensions.
71 
72     /// Modifies the current tensor to have the number and size of
73     /// dimensions as described in the \c vector \c d .  The total number
74     /// of elements must be the same before and after, and the current
75     /// tensor must be contiguous.
reshape_inplace(const std::vector<long> & d)76     void BaseTensor::reshape_inplace(const std::vector<long>& d) {
77         reshape_inplace(d.size(), &d[0]);
78     }
79 
80     /// Reshape the current tensor to be the same size and 1-d. It must be contiguous.
flat_inplace()81     void BaseTensor::flat_inplace() {
82         TENSOR_ASSERT(iscontiguous(),"not contiguous",0,this);
83         long d[] = {_size};
84         set_dims_and_size(1,d);
85     }
86 
87     /// Split dimension i in two ... the product of the new dimensions must match the old.
splitdim_inplace(long i,long dimi0,long dimi1)88     void BaseTensor::splitdim_inplace(long i, long dimi0, long dimi1) {
89         if (i < 0) i += _ndim;
90         TENSOR_ASSERT(i>=0 && i<_ndim, "invalid dimension", i, this);
91         TENSOR_ASSERT(dimi0*dimi1 == _dim[i], "before & after sizes do not match",
92                       _dim[i], this);
93         TENSOR_ASSERT(_ndim+1 <= TENSOR_MAXDIM, "resulting tensor has too many dimensions",
94                       _ndim+1, this);
95         for (long j=_ndim-1; j>i; --j) {
96             _dim[j+1] = _dim[j];
97             _stride[j+1] = _stride[j];
98         }
99         _dim[i+1] = dimi1;
100         _stride[i+1] = _stride[i];
101         _dim[i] = dimi0;
102         _stride[i] *= dimi1;
103         ++_ndim;
104     }
105 
106     /// Fuse the contiguous dimensions i and i+1.
fusedim_inplace(long i)107     void BaseTensor::fusedim_inplace(long i) { // fuse i,i+1 --> i
108         if (i < 0) i += _ndim;
109         TENSOR_ASSERT(i>=0 && i<(_ndim-1) && _ndim>1,"invalid dimension",i,this);
110         TENSOR_ASSERT(_stride[i] == _dim[i+1]*_stride[i+1],"dimensions are not contiguous",
111                       i, this);
112         _dim[i] *= _dim[i+1];
113         _stride[i] = _stride[i+1];
114         for (long j=i+1; j<=_ndim-1; ++j) {
115             _dim[j] = _dim[j+1];
116             _stride[j] = _stride[j+1];
117         }
118         _ndim--;
119         _dim[_ndim] = 1;		// So can iterate over missing dimensions
120         _stride[_ndim] = 0;
121     }
122 
123     /// Swap the dimensions i and j.
swapdim_inplace(long i,long j)124     void BaseTensor::swapdim_inplace(long i, long j) {
125         if (i < 0) i += _ndim;
126         if (j < 0) j += _ndim;
127         TENSOR_ASSERT(i>=0 && i<_ndim,"invalid dimension i",i,this);
128         TENSOR_ASSERT(j>=0 && j<_ndim,"invalid dimension j",j,this);
129         std::swap<long>(_dim[i],_dim[j]);
130         std::swap<long>(_stride[i],_stride[j]);
131     }
132 
133     /// Cyclic shift by nshift places of the inclusive range of dimensions [start,....,end]
cycledim_inplace(long nshift,long start,long end)134     void BaseTensor::cycledim_inplace(long nshift, long start, long end) {
135         long ndshift, dimtmp[TENSOR_MAXDIM], _stridetmp[TENSOR_MAXDIM];
136         if (start < 0) start += _ndim; // Same convention for -ve as in Slice
137         if (end < 0) end += _ndim;
138         TENSOR_ASSERT(start>=0 && start<_ndim,"invalid start dimension",start,this);
139         TENSOR_ASSERT(end>=0 && end>=start,"invalid end dimension",end,this);
140 
141         ndshift = end - start + 1;
142         for (long i=start; i<=end; ++i) {
143             dimtmp[i] = _dim[i];
144             _stridetmp[i] = _stride[i];
145         }
146         for (long i=end; i>=start; --i) {
147             long j = i + nshift;
148             while (j > end) j -= ndshift;
149             while (j < start) j += ndshift;
150             _dim[j] = dimtmp[i];
151             _stride[j] = _stridetmp[i];
152         }
153     }
154 
155     /// General permuation of the dimensions
mapdim_inplace(const std::vector<long> & map)156     void BaseTensor::mapdim_inplace(const std::vector<long>& map) {
157         TENSOR_ASSERT(_ndim == (int) map.size(),"map[] must include all dimensions",map.size(),this);
158         long tmpd[TENSOR_MAXDIM], tmps[TENSOR_MAXDIM];
159         for (long i=0; i<_ndim; ++i) {
160             tmpd[map[i]] = _dim[i];
161             tmps[map[i]] = _stride[i];
162         }
163         for (long i=0; i<_ndim; ++i) {
164             _dim[i] = tmpd[i];
165             _stride[i] = tmps[i];
166         }
167     }
168 }
169