tensor_macros.h

Go to the documentation of this file.

/*
  This file is part of MADNESS.

  Copyright (C) 2007,2010 Oak Ridge National Laboratory

  This program is free software; you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation; either version 2 of the License, or
  (at your option) any later version.

  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  GNU General Public License for more details.

  You should have received a copy of the GNU General Public License
  along with this program; if not, write to the Free Software
  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA

  For more information please contact:

  Robert J. Harrison
  Oak Ridge National Laboratory
  One Bethel Valley Road
  P.O. Box 2008, MS-6367

  email: harrisonrj@ornl.gov
  tel:   865-241-3937
  fax:   865-572-0680


  $Id$
*/


#ifndef MADNESS_TENSOR_TENSOR_MACROS_H__INCLUDED
#define MADNESS_TENSOR_TENSOR_MACROS_H__INCLUDED

// don't change this without changing the iterator macros
#define TENSOR_MAXDIM 6


#define IND1 _i
#define IND2 _i,_j
#define IND3 _i,_j,_k
#define IND4 _i,_j,_k,_l
#define IND5 _i,_j,_k,_l,_m
#define IND6 _i,_j,_k,_l,_m,_n
#define IND  IND6


#define ITERATOR1(t,exp) do { \
        long __xd0=t.dim(0),_index=0;                                   \
        for (long _i=0; _i<__xd0; ++_i) {exp;_index++;} } while (0)

#define ITERATOR2(t,exp) do { \
        long __xd0=t.dim(0), __xd1=t.dim(1), _index=0;  \
for (long _i=0; _i<__xd0; ++_i) { \
  for (long _j=0; _j<__xd1; ++_j) {exp;_index++;} } } while (0)

#define ITERATOR3(t,exp) do { \
        long __xd0=t.dim(0), __xd1=t.dim(1), __xd2=t.dim(2), _index=0;  \
for (long _i=0; _i<__xd0; ++_i) { \
  for (long _j=0; _j<__xd1; ++_j) { \
    for (long _k=0; _k<__xd2; ++_k) {exp;_index++;} } } } while (0)

#define ITERATOR4(t,exp) do { \
        long __xd0=t.dim(0), __xd1=t.dim(1), __xd2=t.dim(2),        \
            __xd3=t.dim(3), _index=0;                               \
for (long _i=0; _i<__xd0; ++_i) { \
  for (long _j=0; _j<__xd1; ++_j) { \
    for (long _k=0; _k<__xd2; ++_k) { \
      for (long _l=0; _l<__xd3; ++_l) {exp;_index++;} } } } } while (0)

#define ITERATOR5(t,exp) do { \
        long __xd0=t.dim(0), __xd1=t.dim(1), __xd2=t.dim(2),      \
            __xd3=t.dim(3), __xd4=t.dim(4), _index=0;             \
for (long _i=0; _i<__xd0; ++_i) { \
  for (long _j=0; _j<__xd1; ++_j) { \
    for (long _k=0; _k<__xd2; ++_k) { \
      for (long _l=0; _l<__xd3; ++_l) { \
        for (long _m=0; _m<__xd4; ++_m) {exp;_index++;} } } } } } while (0)

#define ITERATOR6(t,exp) do { \
        long __xd0=t.dim(0), __xd1=t.dim(1), __xd2=t.dim(2),            \
            __xd3=t.dim(3), __xd4=t.dim(4), __xd5=t.dim(5), _index=0;;  \
for (long _i=0; _i<__xd0; ++_i) { \
  for (long _j=0; _j<__xd1; ++_j) { \
    for (long _k=0; _k<__xd2; ++_k) { \
      for (long _l=0; _l<__xd3; ++_l) { \
        for (long _m=0; _m<__xd4; ++_m) { \
          for (long _n=0; _n<__xd5; ++_n) {exp;_index++;} } } } } } } while(0)

#define ITERATOR(t,exp) do { \
  long _j=0, _k=0, _l=0, _m=0, _n=0; \
  if (t.ndim() == 1) {ITERATOR1(t,exp);} \
  else if (t.ndim() == 2) {ITERATOR2(t,exp);} \
  else if (t.ndim() == 3) {ITERATOR3(t,exp);} \
  else if (t.ndim() == 4) {ITERATOR4(t,exp);} \
  else if (t.ndim() == 5) {ITERATOR5(t,exp);} \
  else if (t.ndim() == 6) {ITERATOR6(t,exp);} \
  else {TENSOR_ASSERT(t.ndim() <= 6,"ndim confused?",t.ndim(),&t);} \
 } while(0)

// Inside iterator access pointer to current element as _p0 (pointer to
// argument number 0).  _i, _j, _k, ..., also defined
#define UNARYITERATOR1(X,x,exp) do { \
        long __xd0=x.dim(0);         \
        long __xs0=x.stride(0);      \
X* restrict _p0=x.ptr(); \
for (long _i=0; _i<__xd0; ++_i,_p0+=__xs0) { \
  exp; \
} } while(0)

#define UNARYITERATOR2(X,x,exp) do { \
        long __xd0=x.dim(0), __xd1=x.dim(1);   \
        long __xs0=x.stride(0), __xs1=x.stride(1);      \
X* restrict __xp0=x.ptr(); \
for (long _i=0; _i<__xd0; ++_i,__xp0+=__xs0) { \
  X* restrict _p0=__xp0; \
  for (long _j=0; _j<__xd1; ++_j, _p0+=__xs1) { \
    exp; \
  } } } while(0)

#define UNARYITERATOR3(X,x,exp) do { \
        long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2);        \
        long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2);   \
X* restrict __xp0=x.ptr(); \
for (long _i=0; _i<__xd0; ++_i,__xp0+=__xs0) { \
  X* restrict __xp1=__xp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1) { \
    X* restrict _p0=__xp1; \
    for (long _k=0; _k<__xd2; ++_k, _p0+=__xs2) { \
       exp; \
    } } } } while(0)

#define UNARYITERATOR4(X,x,exp) do { \
        long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),    \
            __xd3=x.dim(3);                                         \
        long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),   \
            __xs3=x.stride(3);                                          \
X* restrict __xp0=x.ptr(); \
for (long _i=0; _i<__xd0; ++_i,__xp0+=__xs0) { \
  X* restrict __xp1=__xp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1) { \
    X* restrict __xp2=__xp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2) { \
       X* restrict _p0=__xp2; \
       for (long _l=0; _l<__xd3; ++_l, _p0+=__xs3) { \
          exp; \
       } } } } } while(0)

#define UNARYITERATOR5(X,x,exp) do { \
        long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),    \
            __xd3=x.dim(3), __xd4=x.dim(4);                         \
        long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),   \
            __xs3=x.stride(3), __xs4=x.stride(4);                       \
X* restrict __xp0=x.ptr(); \
for (long _i=0; _i<__xd0; ++_i,__xp0+=__xs0) { \
  X* restrict __xp1=__xp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1) { \
    X* restrict __xp2=__xp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2) { \
       X* restrict __xp3=__xp2; \
       for (long _l=0; _l<__xd3; ++_l, __xp3+=__xs3) { \
          X* restrict _p0 =__xp3; \
          for (long _m=0; _m<__xd4; ++_m, _p0+=__xs4) { \
            exp; \
          } } } } } } while(0)

#define UNARYITERATOR6(X,x,exp) do { \
        long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),    \
            __xd3=x.dim(3), __xd4=x.dim(4), __xd5=x.dim(5);         \
        long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),   \
            __xs3=x.stride(3), __xs4=x.stride(4), __xs5=x.stride(5);    \
X* restrict __xp0=x.ptr(); \
for (long _i=0; _i<__xd0; ++_i,__xp0+=__xs0) { \
  X* restrict __xp1=__xp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1) { \
    X* restrict __xp2=__xp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2) { \
       X* restrict __xp3=__xp2; \
       for (long _l=0; _l<__xd3; ++_l, __xp3+=__xs3) { \
          X* restrict __xp4=__xp3; \
          for (long _m=0; _m<__xd4; ++_m, __xp4+=__xs4) { \
            X* restrict _p0=__xp4; \
            for (long _n=0; _n<__xd5; ++_n, _p0+=__xs5) { \
              exp; \
          } } } } } } } while(0)

#define UNARYITERATOR(X,x,exp) do { \
  long _j=0, _k=0, _l=0, _m=0, _n=0; \
  if (x.ndim() == 1) UNARYITERATOR1(X,x,exp); \
  else if (x.ndim() == 2) UNARYITERATOR2(X,x,exp); \
  else if (x.ndim() == 3) UNARYITERATOR3(X,x,exp); \
  else if (x.ndim() == 4) UNARYITERATOR4(X,x,exp); \
  else if (x.ndim() == 5) UNARYITERATOR5(X,x,exp); \
  else if (x.ndim() == 6) UNARYITERATOR6(X,x,exp); \
  else {TENSOR_ASSERT(x.ndim() <= 6,"ndim confused?",x.ndim(),&x);} } while(0)

// Inside iterator access pointers to current elements as _p0 & _p1
// _i, _j, _k, ... also defined
#define BINARYITERATOR1(X,x,Y,y,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
 long __xd0=x.dim(0);                                                   \
 long __xs0=x.stride(0);                                                \
 long __ys0=y.stride(0);                                                \
X* restrict _p0=x.ptr(); \
Y* restrict _p1=y.ptr(); \
for (long _i=0; _i<__xd0; ++_i, _p0+=__xs0, _p1+=__ys0) { \
   exp; \
} } while(0)

#define BINARYITERATOR2(X,x,Y,y,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1);                                   \
 long __xs0=x.stride(0), __xs1=x.stride(1);                             \
 long __ys0=y.stride(0), __ys1=y.stride(1);                             \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0) { \
  X* restrict _p0=__xp0; \
  Y* restrict _p1=__yp0; \
  for (long _j=0; _j<__xd1; ++_j, _p0+=__xs1, _p1+=__ys1) { \
     exp; \
  } } } while(0)

#define BINARYITERATOR3(X,x,Y,y,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2);                   \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2);          \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2);          \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1) { \
    X* restrict _p0=__xp1; \
    Y* restrict _p1=__yp1; \
    for (long _k=0; _k<__xd2; ++_k, _p0+=__xs2, _p1+=__ys2) { \
       exp; \
    } } } } while(0)

#define BINARYITERATOR4(X,x,Y,y,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),                   \
     __xd3=x.dim(3);                                                    \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),          \
     __xs3=x.stride(3);                                                 \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2),          \
     __ys3=y.stride(3);                                                 \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1) { \
    X* restrict __xp2=__xp1; \
    Y* restrict __yp2=__yp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2, __yp2+=__ys2) { \
      X* restrict _p0=__xp2; \
      Y* restrict _p1=__yp2; \
      for (long _l=0; _l<__xd3; ++_l, _p0+=__xs3, _p1+=__ys3) { \
         exp; \
      } } } } } while(0)

#define BINARYITERATOR5(X,x,Y,y,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),                   \
     __xd3=x.dim(3), __xd4=x.dim(4);                                    \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),          \
     __xs3=x.stride(3), __xs4=x.stride(4);                              \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2),          \
     __ys3=y.stride(3), __ys4=y.stride(4);                              \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1) { \
    X* restrict __xp2=__xp1; \
    Y* restrict __yp2=__yp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2, __yp2+=__ys2) { \
      X* restrict __xp3=__xp2; \
      Y* restrict __yp3=__yp2; \
      for (long _l=0; _l<__xd3; ++_l, __xp3+=__xs3, __yp3+=__ys3) { \
        X* restrict _p0=__xp3; \
        Y* restrict _p1=__yp3; \
        for (long _m=0; _m<__xd4; ++_m, _p0+=__xs4, _p1+=__ys4) { \
           exp; \
        } } } } } } while(0)

#define BINARYITERATOR6(X,x,Y,y,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),                   \
     __xd3=x.dim(3), __xd4=x.dim(4), __xd5=x.dim(5);                    \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),          \
     __xs3=x.stride(3), __xs4=x.stride(4), __xs5=x.stride(5);           \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2),          \
     __ys3=y.stride(3), __ys4=y.stride(4), __ys5=y.stride(5);           \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1) { \
    X* restrict __xp2=__xp1; \
    Y* restrict __yp2=__yp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2, __yp2+=__ys2) { \
      X* restrict __xp3=__xp2; \
      Y* restrict __yp3=__yp2; \
      for (long _l=0; _l<__xd3; ++_l, __xp3+=__xs3, __yp3+=__ys3) { \
        X* restrict __xp4=__xp3; \
        Y* restrict __yp4=__yp3; \
        for (long _m=0; _m<__xd4; ++_m, __xp4+=__xs4, __yp4+=__ys4) { \
          X* restrict _p0=__xp4; \
          Y* restrict _p1=__yp4; \
          for (long _n=0; _n<__xd5; ++_n, _p0+=__xs5, _p1+=__ys5) { \
             exp; \
          } } } } } } } while(0)

#define BINARYITERATOR(X,x,Y,y,exp) do { \
  long _j=0, _k=0, _l=0, _m=0, _n=0; \
  if (x.ndim() == 1) BINARYITERATOR1(X,x,Y,y,exp); \
  else if (x.ndim() == 2) BINARYITERATOR2(X,x,Y,y,exp); \
  else if (x.ndim() == 3) BINARYITERATOR3(X,x,Y,y,exp); \
  else if (x.ndim() == 4) BINARYITERATOR4(X,x,Y,y,exp); \
  else if (x.ndim() == 5) BINARYITERATOR5(X,x,Y,y,exp); \
  else if (x.ndim() == 6) BINARYITERATOR6(X,x,Y,y,exp); \
  else {TENSOR_ASSERT(x.ndim() <= 6,"ndim confused?",x.ndim(),&x);} \
} while(0)

// Inside iterator access pointers to current elements as _p0, _p1, _p2
// _i, _j, _k, ... also defined
#define TERNARYITERATOR1(X,x,Y,y,Z,z,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
TENSOR_ASSERT(x.conforms(z),"first and third tensors do not conform",0,&x); \
 long __xd0=x.dim(0);                                                   \
 long __xs0=x.stride(0);                                                \
 long __ys0=y.stride(0);                                                \
 long __zs0=z.stride(0);                                                \
X* restrict _p0=x.ptr(); \
Y* restrict _p1=y.ptr(); \
Z* restrict _p2=z.ptr(); \
for (long _i=0; _i<__xd0; ++_i, _p0+=__xs0, _p1+=__ys0, _p2+=__zs0) { \
  exp; \
} } while(0)

#define TERNARYITERATOR2(X,x,Y,y,Z,z,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
TENSOR_ASSERT(x.conforms(z),"first and third tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1);                                   \
 long __xs0=x.stride(0), __xs1=x.stride(1);                             \
 long __ys0=y.stride(0), __ys1=y.stride(1);                             \
 long __zs0=z.stride(0), __zs1=z.stride(1);                             \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
Z* restrict __zp0=z.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0, __zp0+=__zs0) { \
  X* restrict _p0=__xp0; \
  Y* restrict _p1=__yp0; \
  Z* restrict _p2=__zp0; \
  for (long _j=0; _j<__xd1; ++_j, _p0+=__xs1, _p1+=__ys1, _p2+=__zs1) { \
    exp; \
  } } } while(0)

#define TERNARYITERATOR3(X,x,Y,y,Z,z,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
TENSOR_ASSERT(x.conforms(z),"first and third tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2);                   \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2);          \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2);          \
 long __zs0=z.stride(0), __zs1=z.stride(1), __zs2=z.stride(2);          \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
Z* restrict __zp0=z.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0, __zp0+=__zs0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  Z* restrict __zp1=__zp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1, __zp1+=__zs1) { \
    X* restrict _p0=__xp1; \
    Y* restrict _p1=__yp1; \
    Z* restrict _p2=__zp1; \
    for (long _k=0; _k<__xd2; ++_k, _p0+=__xs2, _p1+=__ys2, _p2+=__zs2) { \
      exp; \
    } } } } while(0)

#define TERNARYITERATOR4(X,x,Y,y,Z,z,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
TENSOR_ASSERT(x.conforms(z),"first and third tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),                   \
     __xd3=x.dim(3);                                                    \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),          \
     __xs3=x.stride(3);                                                 \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2),          \
     __ys3=y.stride(3);                                                 \
 long __zs0=z.stride(0), __zs1=z.stride(1), __zs2=z.stride(2),          \
     __zs3=z.stride(3);                                                 \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
Z* restrict __zp0=z.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0, __zp0+=__zs0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  Z* restrict __zp1=__zp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1, __zp1+=__zs1) { \
    X* restrict __xp2=__xp1; \
    Y* restrict __yp2=__yp1; \
    Z* restrict __zp2=__zp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2, __yp2+=__ys2, __zp2+=__zs2) { \
      X* restrict _p0=__xp2; \
      Y* restrict _p1=__yp2; \
      Z* restrict _p2=__zp2; \
      for (long _l=0; _l<__xd3; ++_l, _p0+=__xs3, _p1+=__ys3, _p2+=__zs3) { \
        exp; \
      } } } } } while(0)

#define TERNARYITERATOR5(X,x,Y,y,Z,z,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
TENSOR_ASSERT(x.conforms(z),"first and third tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),                   \
     __xd3=x.dim(3), __xd4=x.dim(4);                                    \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),          \
     __xs3=x.stride(3), __xs4=x.stride(4);                              \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2),          \
     __ys3=y.stride(3), __ys4=y.stride(4);                              \
 long __zs0=z.stride(0), __zs1=z.stride(1), __zs2=z.stride(2),          \
     __zs3=z.stride(3), __zs4=z.stride(4);                              \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
Z* restrict __zp0=z.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0, __zp0+=__zs0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  Z* restrict __zp1=__zp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1, __zp1+=__zs1) { \
    X* restrict __xp2=__xp1; \
    Y* restrict __yp2=__yp1; \
    Z* restrict __zp2=__zp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2, __yp2+=__ys2, __zp2+=__zs2) { \
      X* restrict __xp3=__xp2; \
      Y* restrict __yp3=__yp2; \
      Z* restrict __zp3=__zp2; \
      for (long _l=0; _l<__xd3; ++_l, __xp3+=__xs3, __yp3+=__ys3, __zp3+=__zs3) { \
        X* restrict _p0=__xp3; \
        Y* restrict _p1=__yp3; \
        Z* restrict _p2=__zp3; \
        for (long _m=0; _m<__xd4; ++_m, _p0+=__xs4, _p1+=__ys4, _p2+=__zs4) { \
          exp; \
        } } } } } } while(0)

#define TERNARYITERATOR6(X,x,Y,y,Z,z,exp) do { \
TENSOR_ASSERT(x.conforms(y),"first and second tensors do not conform",0,&x); \
TENSOR_ASSERT(x.conforms(z),"first and third tensors do not conform",0,&x); \
 long __xd0=x.dim(0), __xd1=x.dim(1), __xd2=x.dim(2),                   \
     __xd3=x.dim(3), __xd4=x.dim(4), __xd5=x.dim(5);                    \
 long __xs0=x.stride(0), __xs1=x.stride(1), __xs2=x.stride(2),          \
     __xs3=x.stride(3), __xs4=x.stride(4), __xs5=x.stride(5);           \
 long __ys0=y.stride(0), __ys1=y.stride(1), __ys2=y.stride(2),          \
     __ys3=y.stride(3), __ys4=y.stride(4), __ys5=y.stride(5);           \
 long __zs0=z.stride(0), __zs1=z.stride(1), __zs2=z.stride(2),          \
     __zs3=z.stride(3), __zs4=z.stride(4), __zs5=z.stride(5);           \
X* restrict __xp0=x.ptr(); \
Y* restrict __yp0=y.ptr(); \
Z* restrict __zp0=z.ptr(); \
for (long _i=0; _i<__xd0; ++_i, __xp0+=__xs0, __yp0+=__ys0, __zp0+=__zs0) { \
  X* restrict __xp1=__xp0; \
  Y* restrict __yp1=__yp0; \
  Z* restrict __zp1=__zp0; \
  for (long _j=0; _j<__xd1; ++_j, __xp1+=__xs1, __yp1+=__ys1, __zp1+=__zs1) { \
    X* restrict __xp2=__xp1; \
    Y* restrict __yp2=__yp1; \
    Z* restrict __zp2=__zp1; \
    for (long _k=0; _k<__xd2; ++_k, __xp2+=__xs2, __yp2+=__ys2, __zp2+=__zs2) { \
      X* restrict __xp3=__xp2; \
      Y* restrict __yp3=__yp2; \
      Z* restrict __zp3=__zp2; \
      for (long _l=0; _l<__xd3; ++_l, __xp3+=__xs3, __yp3+=__ys3, __zp3+=__zs3) { \
        X* restrict __xp4=__xp3; \
        Y* restrict __yp4=__yp3; \
        Z* restrict __zp4=__zp3; \
        for (long _m=0; _m<__xd4; ++_m, __xp4+=__xs4, __yp4+=__ys4, __zp4+=__zs4) { \
          X* restrict _p0=__xp4; \
          Y* restrict _p1=__yp4; \
          Z* restrict _p2=__zp4; \
          for (long _n=0; _n<__xd5; ++_n, _p0+=__xs5, _p1+=__ys5, _p2+=__zs5) { \
            exp; \
          } } } } } } } while(0)

#define TERNARYITERATOR(X,x,Y,y,Z,z,exp) do { \
  long _j=0, _k=0, _l=0, _m=0, _n=0; \
  if (x.ndim() == 1) TERNARYITERATOR1(X,x,Y,y,Z,z,exp); \
  else if (x.ndim() == 2) TERNARYITERATOR2(X,x,Y,y,Z,z,exp); \
  else if (x.ndim() == 3) TERNARYITERATOR3(X,x,Y,y,Z,z,exp); \
  else if (x.ndim() == 4) TERNARYITERATOR4(X,x,Y,y,Z,z,exp); \
  else if (x.ndim() == 5) TERNARYITERATOR5(X,x,Y,y,Z,z,exp); \
  else if (x.ndim() == 6) TERNARYITERATOR6(X,x,Y,y,Z,z,exp); \
  else {TENSOR_ASSERT(x.ndim() <= 6,"ndim confused?",x.ndim(),&x);} \
} while(0)

#define UNARY_OPTIMIZED_ITERATOR(X,x,exp) do { \
  if (x.iscontiguous()) { \
    X* restrict _p0 = x.ptr(); \
    for (long _j=0; _j<x.size(); ++_j,++_p0) {exp;} \
  } \
  else { \
    for (TensorIterator<REMCONST(X)> iter=x.unary_iterator(1); iter._p0; ++iter) { \
      long _dimj = iter.dimj; \
      X* restrict _p0 = iter._p0; \
      long _s0 = iter._s0; \
      for (long _j=0; _j<_dimj; ++_j, _p0+=_s0) { \
        exp; \
      } \
    } \
  } \
} while(0)


// Can optimize these by moving definition of stride out of loop.

#define UNARY_UNOPTIMIZED_ITERATOR(X,x,exp) do { \
    for (TensorIterator<REMCONST(X)> iter=x.unary_iterator(1,false,false); iter._p0; ++iter) { \
    long _dimj = iter.dimj; \
    X* restrict _p0 = iter._p0; \
    long _s0 = iter._s0; \
    for (long _j=0; _j<_dimj; ++_j, _p0+=_s0) { \
      exp; \
    } \
  } } while(0)

// Use this inside another iterator macro ... pointers are _q* instead of _p*
// Iterator is iter2.

#define UNARY_UNOPTIMIZED_ITERATOR_NESTED(X,x,exp) do { \
    for (TensorIterator<REMCONST(X)> iter2=x.unary_iterator(1,false,false); iter2._p0; ++iter2) { \
    long _dimj2 = iter2.dimj; \
    X* restrict _q0 = iter2._p0; \
    long _s20 = iter2._s0; \
    for (long _j2=0; _j2<_dimj2; ++_j2, _q0+=_s20) { \
      exp; \
    } \
  } } while(0)

#define BINARY_OPTIMIZED_ITERATOR(X,x,Y,y,exp) do { \
  if (x.iscontiguous() && y.iscontiguous() && x.size()==y.size()) { \
    X* restrict _p0 = x.ptr(); \
    Y* restrict _p1 = y.ptr(); \
    for (long _j=0; _j<x.size(); ++_j,++_p0,++_p1) {exp;} \
  } \
  else { \
    for (TensorIterator<REMCONST(X),REMCONST(Y)> iter=x.binary_iterator(y,1); iter._p0; ++iter) { \
        long _dimj = iter.dimj; \
        X* restrict _p0 = iter._p0; \
        Y* restrict _p1 = iter._p1; \
        long _s0 = iter._s0; \
        long _s1 = iter._s1; \
        for (long _j=0; _j<_dimj; ++_j, _p0+=_s0, _p1+=_s1) { \
          exp; \
        } \
  } } } while(0)

#define TERNARY_OPTIMIZED_ITERATOR(X,x,Y,y,Z,z,exp) do { \
  if (x.iscontiguous() && y.iscontiguous() && z.iscontiguous() && x.size()==y.size() && x.size()==z.size()) { \
    X* restrict _p0 = x.ptr(); \
    Y* restrict _p1 = y.ptr(); \
    Z* restrict _p2 = z.ptr(); \
    for (long _j=0; _j<x.size(); ++_j,++_p0,++_p1,++_p2) {exp;} \
  } \
  else { \
    for (TensorIterator<REMCONST(X),REMCONST(Y),REMCONST(Z)> iter=x.ternary_iterator(y,z,1); iter._p0; ++iter) { \
        long _dimj = iter.dimj; \
        X* restrict _p0 = iter._p0; \
        Y* restrict _p1 = iter._p1; \
        Z* restrict _p2 = iter._p2; \
        long _s0 = iter._s0; \
        long _s1 = iter._s1; \
        long _s2 = iter._s2; \
        for (long _j=0; _j<_dimj; ++_j, _p0+=_s0, _p1+=_s1, _p2+=_s2) { \
          exp; \
        } \
  } } } while(0)

#endif // MADNESS_TENSOR_TENSOR_MACROS_H__INCLUDED