   }
 };
 
+
+
+
+template<class DeviceType, class LeftViewType, class RightViewType, class OutputViewType>
+struct ContractFieldFieldTensorFunctor1D {
+  typedef DeviceType device_type;
+  LeftViewType _leftFields;
+  RightViewType _rightFields;
+  OutputViewType _outputFields;
+  int _numLeftFields;
+  int _numRightFields;
+  int _numPoints;
+  int _dim1Tensor;
+  int _dim2Tensor;
+  int _numCells;
+
+  ContractFieldFieldTensorFunctor1D(LeftViewType leftFields,
+				  RightViewType rightFields,
+				  OutputViewType outputFields,
+				  int numLeftFields,
+				  int numRightFields,
+				  int numPoints,
+				  int dim1Tensor,
+				  int dim2Tensor,
+				  int numCells) :
+    _leftFields(leftFields),
+    _rightFields(rightFields),
+    _outputFields(outputFields),
+    _numLeftFields(numLeftFields),
+    _numRightFields(numRightFields),
+    _numPoints(numPoints),
+    _dim1Tensor(dim1Tensor),
+    _dim2Tensor(dim2Tensor),
+    _numCells(numCells)
+  {
+    // Nothing to do
+  }
+
+  KOKKOS_INLINE_FUNCTION
+  void operator()(const unsigned int elementIndex) const {
+
+    for (int lbf = 0; lbf < _numLeftFields; lbf++) {
+      for (int rbf = 0; rbf < _numRightFields; rbf++) {
+        double tmpVal = 0;
+        for (int qp = 0; qp < _numPoints; qp++) {
+          for (int iTens1 = 0; iTens1 < _dim1Tensor; iTens1++) {
+            for (int iTens2 = 0; iTens2 < _dim2Tensor; iTens2++) {
+              tmpVal +=
+	      _leftFields(lbf*_numPoints*_dim1Tensor*_dim2Tensor*_numCells +
+	      qp*_dim1Tensor*_dim2Tensor*_numCells +
+	      iTens1*_dim2Tensor*_numCells + iTens2*_numCells + elementIndex)
+	      *_rightFields(rbf*_numPoints*_dim1Tensor*_dim2Tensor*_numCells +
+	      qp*_dim1Tensor*_dim2Tensor*_numCells +
+	      iTens1*_dim2Tensor*_numCells + iTens2*_numCells + elementIndex);
+            } // D2-loop
+          } // D1-loop
+        } // P-loop
+        _outputFields(lbf*_numRightFields*_numCells +
+	rbf*_numCells + elementIndex) = tmpVal;
+      } // R-loop
+    } // L-loop
+  }
+};
+
+
+
+
+
 // Serial contractFieldFieldTensor.  Contracts FieldContainers of doubles.
 void contractFieldFieldTensorSerial(FieldContainer<double> &  outputFields,
 			      const FieldContainer<double> &              leftFields,

  double elapsedTime_kokkos_noCopy = 0;

  for(int i = 0; i < 5; ++i){
    contractFieldFieldTensorKokkos<Kokkos::OpenMP, omp_input_view_t, omp_output_view_t, omp_input_host_t, omp_output_host_t> 
      (omp_hostOut, omp_hostLeft, omp_hostRight, omp_kokkosOut, omp_kokkosLeft,
       omp_kokkosRight, &elapsedTime_kokkos_noCopy);
  }

  // Copy out from kokkos output view (NOT timing this)
  for (int cl = 0; cl < c; ++cl) {
    for (int lbf = 0; lbf < l; ++lbf) {
      for (int rbf = 0; rbf < r; ++rbf) {
        out1_c_l_r(cl, lbf, rbf) = omp_hostOut(cl, lbf, rbf);
      }
    }
  }

  // device.

  Kokkos::initialize();
  
  // Kokkos Cuda views
  typedef Kokkos::View<double*, Kokkos::Cuda> cuda_input_view_t;
  typedef Kokkos::View<double *, Kokkos::Cuda> cuda_input_view_t;
  typedef Kokkos::View<double *, Kokkos::Cuda> cuda_output_view_t;
  typedef typename cuda_input_view_t::HostMirror cuda_input_host_t;
  typedef typename cuda_output_view_t::HostMirror cuda_output_host_t;

  // Kokkos OpenMP views
  typedef Kokkos::View<double *, Kokkos::OpenMP> omp_input_view_t;
  typedef Kokkos::View<double *, Kokkos::OpenMP> omp_output_view_t;
  typedef typename omp_input_view_t::HostMirror omp_input_host_t;
  typedef typename omp_output_view_t::HostMirror omp_output_host_t;

  // Make equivalent Kokkos views
  cuda_input_view_t cuda_kokkosLeft("left_input", c * l * p * d1 * d2);
  cuda_input_view_t cuda_kokkosRight("right_input", c * r * p * d1 * d2);
  cuda_output_view_t cuda_kokkosOut("output", c * l * r);

  omp_input_view_t omp_kokkosLeft("left_input", c * l * p * d1 * d2);
  omp_input_view_t omp_kokkosRight("right_input",  c * r * p * d1 * d2);
  omp_output_view_t omp_kokkosOut("output", c * l * r);

  // And their host mirrors
  cuda_input_host_t cuda_hostLeft = Kokkos::create_mirror_view(cuda_kokkosLeft);
  cuda_input_host_t cuda_hostRight = Kokkos::create_mirror_view(cuda_kokkosRight);
  cuda_output_host_t cuda_hostOut = Kokkos::create_mirror_view(cuda_kokkosOut);

  omp_input_host_t omp_hostLeft = Kokkos::create_mirror_view(omp_kokkosLeft);
  omp_input_host_t omp_hostRight = Kokkos::create_mirror_view(omp_kokkosRight);
  omp_output_host_t omp_hostOut = Kokkos::create_mirror_view(omp_kokkosOut);
  
  // Copy into Kokkos host views.
  // Need to change this so that its 1-D and cl has stride 1
  for (int cl = 0; cl < c; ++cl) {
    for (int lbf = 0; lbf < l; ++lbf) {
      for (int qp = 0; qp < p; ++qp) {
        for (int iTens1 = 0; iTens1 < d1; ++iTens1) {
          for (int iTens2 = 0; iTens2 < d2; ++iTens2) {
            cuda_hostLeft(lbf*p*d1*d2*c + qp*d1*d2*c + iTens1*d2*c + iTens2*c +
	    cl) = in_c_l_p_d_d(cl, lbf, qp, iTens1, iTens2);
            omp_hostLeft(lbf*p*d1*d2*c + qp*d1*d2*c + iTens1*d2*c + iTens2*c +
	    cl) = in_c_l_p_d_d(cl, lbf, qp, iTens1, iTens2);
          }
        }
      }
    }
    for (int rbf = 0; rbf < r; ++rbf) {
      for (int qp = 0; qp < p; ++qp) {
        for (int iTens1 = 0; iTens1 < d1; ++iTens1) {
          for (int iTens2 = 0; iTens2 < d2; ++iTens2) {
            cuda_hostRight(rbf*p*d1*d2*c + qp*d1*d2*c + iTens1*d2*c + iTens2*c +
	    cl) = in_c_r_p_d_d(cl, rbf, qp, iTens1, iTens2);
            omp_hostRight(rbf*p*d1*d2*c + qp*d1*d2*c + iTens1*d2*c + iTens2*c +
	    cl) = in_c_r_p_d_d(cl, rbf, qp, iTens1, iTens2);
          }
        }
      }
    }
  }








#if 0
  // Kokkos Cuda views
  typedef Kokkos::View<double *****, Kokkos::Cuda> cuda_input_view_t;
  typedef Kokkos::View<double ***, Kokkos::Cuda> cuda_output_view_t;
