CUTLASS 3.1 (#915)

Co-authored-by: Aniket Shivam <ashivam@nvidia.com>
2023-04-14 20:19:34 -07:00
parent 9b8166e3f0
commit d572cc1aab
482 changed files with 37184 additions and 16419 deletions
--- a/tools/util/include/cutlass/util/gett_commandline.hpp
+++ b/tools/util/include/cutlass/util/gett_commandline.hpp
@@ -209,7 +209,7 @@ struct GettCommandLine {
    //

    // Permute the batched modes to promote coalescing
-    // Sort the batched modes by min(ldAl,ldBl) and tie-broken by the size
+    // Sort the batched modes by min(ldAl,ldBl) and in case of a tie by the size
    std::sort(std::begin(bat_mode), std::end(bat_mode), [&](char l1, char l2) {
        return std::tie(std::min(mode_ldA[l1],mode_ldB[l1]),mode_size[l1])
             < std::tie(std::min(mode_ldA[l2],mode_ldB[l2]),mode_size[l2]);
@@ -227,7 +227,7 @@ struct GettCommandLine {
    }

    // Permute the reduction modes to promote coalescing
-    // Sort the reduction modes by min(ldAk,ldBk) and tie-broken by the size
+    // Sort the reduction modes by min(ldAk,ldBk) and in case of a tie by the size
    std::sort(std::begin(red_mode), std::end(red_mode), [&](char k1, char k2) {
        return std::tie(std::min(mode_ldA[k1],mode_ldB[k1]),mode_size[k1])
             < std::tie(std::min(mode_ldA[k2],mode_ldB[k2]),mode_size[k2]);
@@ -243,7 +243,7 @@ struct GettCommandLine {
    }

    // Permute the row modes to promote coalescing
-    // Sort the row modes by min(ldAm,ldCm) and tie-broken by ldAm
+    // Sort the row modes by min(ldAm,ldCm) and in case of a tie by ldAm
    std::sort(std::begin(row_mode), std::end(row_mode), [&](char m1, char m2) {
        return std::tie(std::min(mode_ldA[m1],mode_ldC[m1]),mode_ldA[m1])
             < std::tie(std::min(mode_ldA[m2],mode_ldC[m2]),mode_ldA[m2]);
@@ -259,7 +259,7 @@ struct GettCommandLine {
    }

    // Permute the col modes to promote coalescing
-    // Sort the col modes by min(ldBn,ldCn) and tie-broken by ldBn
+    // Sort the col modes by min(ldBn,ldCn) and in case of a tie by ldBn
    std::sort(std::begin(col_mode), std::end(col_mode), [&](char n1, char n2) {
        return std::tie(std::min(mode_ldB[n1],mode_ldC[n1]),mode_ldB[n1])
             < std::tie(std::min(mode_ldB[n2],mode_ldC[n2]),mode_ldB[n2]);
@@ -362,7 +362,7 @@ struct GettCommandLine {
      "    A command delimited list of symbolic mode and its corresponding extent.\n"
      "    Extents are defaulted to 1 if any are not provided.\n\n"

-      "Example usage: gett.exe --modeC=m,n,l --modeA=m,k,l --modeB=k,n,l --extent=m:4096,n:4096,k:4096\n";
+      "Example usage: gett.exe --modeC=m,n,l --modeA=m,k,l --modeB=k,n,l --extents=m:4096,n:4096,k:4096\n";
  }
 };

--- a/tools/util/include/cutlass/util/reference/device/gemm_complex.h
+++ b/tools/util/include/cutlass/util/reference/device/gemm_complex.h
@@ -67,7 +67,8 @@ template <
  typename LayoutC,
  typename ScalarType,
  typename ComputeType,
-  typename ConvertOp = NumericConverter<ElementC, ScalarType>,
+  typename ElementD = ElementC,
+  typename ConvertOp = NumericConverter<ElementD, ScalarType>,
  typename InnerProductOp = multiply_add<ComputeType>,
  int kMblock = 4,
  int kNblock = 4
@@ -81,7 +82,7 @@ __global__ void GemmComplex(
  ComplexTransform transform_b,
  ScalarType beta,
  TensorRef<ElementC, LayoutC> tensor_c,
-  TensorRef<ElementC, LayoutC> tensor_d,
+  TensorRef<ElementD, LayoutC> tensor_d,
  ComputeType initial_accum,
  int batch_count = 1,
  int64_t batch_stride_A = 0,
@@ -198,7 +199,8 @@ template <
  typename LayoutC,
  typename ScalarType,
  typename ComputeType,
-  typename ConvertOp = NumericConverter<ElementC, ScalarType>,
+  typename ElementD = ElementC,
+  typename ConvertOp = NumericConverter<ElementD, ScalarType>,
  typename InnerProductOp = multiply_add<ComputeType>
 >
 void GemmComplex(
@@ -210,7 +212,7 @@ void GemmComplex(
  ComplexTransform transform_b,
  ScalarType beta,
  TensorRef<ElementC, LayoutC> tensor_c,
-  TensorRef<ElementC, LayoutC> tensor_d,
+  TensorRef<ElementD, LayoutC> tensor_d,
  ComputeType initial_accum,
  int batch_count = 1,
  int64_t batch_stride_A = 0,
@@ -243,6 +245,7 @@ void GemmComplex(
      LayoutC,
      ScalarType,
      ComputeType,
+      ElementD,
      ConvertOp,
      InnerProductOp,
      kMblock,
@@ -285,6 +288,7 @@ void GemmComplex(
      LayoutC,
      ScalarType,
      ComputeType,
+      ElementD,
      ConvertOp,
      InnerProductOp,
      kBigMblock,
@@ -322,7 +326,8 @@ template <
  typename LayoutB,
  typename ElementC,
  typename LayoutC,
-  typename ScalarType
+  typename ScalarType,
+  typename ElementD = ElementC
 >
 void GemmComplex(
  gemm::GemmCoord problem_size,
@@ -333,7 +338,7 @@ void GemmComplex(
  ComplexTransform transform_b,
  ScalarType beta,
  TensorRef<ElementC, LayoutC> tensor_c,
-  TensorRef<ElementC, LayoutC> tensor_d) {
+  TensorRef<ElementD, LayoutC> tensor_d) {

  GemmComplex(problem_size, alpha, tensor_a, transform_a, tensor_b, transform_b, beta, tensor_c, tensor_d, ScalarType(0));
 }
--- a/tools/util/include/cutlass/util/reference/device/tensor_fill.h
+++ b/tools/util/include/cutlass/util/reference/device/tensor_fill.h
@@ -58,6 +58,8 @@
 #include "cutlass/tensor_view.h"
 #include "cutlass/blas3.h"

+#include "cutlass/layout/vector.h"
+
 #include "cutlass/util/reference/device/tensor_foreach.h"
 #include "cutlass/util/distribution.h"

@@ -1646,6 +1648,15 @@ void BlockFillSequential(
  Element v = Element(1),
  Element s = Element(0)) {

+  using Layout = layout::PackedVectorLayout;
+  Layout::TensorCoord size(static_cast<Layout::Index>(capacity)); // -Wconversion
+  Layout layout = Layout::packed(size);
+  TensorView<Element, Layout> view(ptr, layout, size);
+
+  Array<Element, Layout::kRank> c;
+  c[0] = v;
+
+  TensorFillLinear(view, c, s);
 }

 ///////////////////////////////////////////////////////////////////////////////////////////////////
--- a/tools/util/include/cutlass/util/reference/host/convolution.h
+++ b/tools/util/include/cutlass/util/reference/host/convolution.h
@@ -65,7 +65,8 @@ template <
  typename LayoutC,
  typename ElementCompute,
  typename ElementAccumulator = ElementCompute,
-  typename ConvertOp = NumericConverter<ElementC, ElementCompute>,
+  typename ElementD = ElementC,
+  typename ConvertOp = NumericConverter<ElementD, ElementCompute>,
  typename InnerProductOp = multiply_add<ElementAccumulator>
 >
 void Conv2dFprop(
@@ -73,7 +74,7 @@ void Conv2dFprop(
  TensorRef<ElementA, LayoutA> tensor_x,
  TensorRef<ElementB, LayoutB> tensor_w,
  TensorRef<ElementC, LayoutC> tensor_y_in,
-  TensorRef<ElementC, LayoutC> tensor_y_out,
+  TensorRef<ElementD, LayoutC> tensor_y_out,
  ElementCompute alpha,
  ElementCompute beta) {

@@ -142,12 +143,13 @@ template <typename ElementA,
          typename LayoutC,
          typename ElementCompute,
          typename ElementAccumulator = ElementCompute,
-          typename ConvertOp = NumericConverter<ElementC, ElementCompute>,
-          typename InnerProductOp = multiply_add<ElementAccumulator> >
+          typename ElementD = ElementC,
+          typename ConvertOp = NumericConverter<ElementD, ElementCompute>,
+          typename InnerProductOp = multiply_add<ElementAccumulator>>
 void Depsep_Fprop(cutlass::TensorView<ElementA, LayoutA> tensor_A,
                  cutlass::TensorView<ElementB, LayoutB> tensor_B,
                  cutlass::TensorView<ElementC, LayoutC> tensor_C,
-                  cutlass::TensorView<ElementC, LayoutC> tensor_D,
+                  cutlass::TensorView<ElementD, LayoutC> tensor_D,
                  ElementCompute alpha,
                  ElementCompute beta,
                  cutlass::Tensor4DCoord padding = cutlass::Tensor4DCoord(),
@@ -208,7 +210,8 @@ template <
  typename LayoutC,
  typename ElementCompute,
  typename ElementAccumulator = ElementCompute,
-  typename ConvertOp = NumericConverter<ElementC, ElementCompute>,
+  typename ElementD = ElementC,
+  typename ConvertOp = NumericConverter<ElementD, ElementCompute>,
  typename InnerProductOp = multiply_add<ElementAccumulator>
 >
 void Conv2dDgrad(
@@ -216,7 +219,7 @@ void Conv2dDgrad(
  TensorRef<ElementA, LayoutA> tensor_dy,
  TensorRef<ElementB, LayoutB> tensor_w,
  TensorRef<ElementC, LayoutC> tensor_dx_in,
-  TensorRef<ElementC, LayoutC> tensor_dx_out,
+  TensorRef<ElementD, LayoutC> tensor_dx_out,
  ElementCompute alpha,
  ElementCompute beta) {

@@ -309,7 +312,8 @@ template <
  typename LayoutC,
  typename ElementCompute,
  typename ElementAccumulator = ElementCompute,
-  typename ConvertOp = NumericConverter<ElementC, ElementCompute>,
+  typename ElementD = ElementC,
+  typename ConvertOp = NumericConverter<ElementD, ElementCompute>,
  typename InnerProductOp = multiply_add<ElementAccumulator>
 >
 void Conv2dWgrad(
@@ -317,7 +321,7 @@ void Conv2dWgrad(
  TensorRef<ElementA, LayoutA> tensor_dy,
  TensorRef<ElementB, LayoutB> tensor_x,
  TensorRef<ElementC, LayoutC> tensor_dw_in,
-  TensorRef<ElementC, LayoutC> tensor_dw_out,
+  TensorRef<ElementD, LayoutC> tensor_dw_out,
  ElementCompute alpha,
  ElementCompute beta) {
  
@@ -389,7 +393,8 @@ template <
  typename LayoutC,
  typename ElementCompute,
  typename ElementAccumulator = ElementCompute,
-  typename ConvertOp = NumericConverter<ElementC, ElementCompute>,
+  typename ElementD = ElementC,
+  typename ConvertOp = NumericConverter<ElementD, ElementCompute>,
  typename InnerProductOp = multiply_add<ElementAccumulator>
 >
 void Conv2d(
@@ -398,7 +403,7 @@ void Conv2d(
  TensorRef<ElementA, LayoutA> tensor_A,
  TensorRef<ElementB, LayoutB> tensor_B,
  TensorRef<ElementC, LayoutC> tensor_C,
-  TensorRef<ElementC, LayoutC> tensor_D,
+  TensorRef<ElementD, LayoutC> tensor_D,
  ElementCompute alpha,
  ElementCompute beta) {

@@ -409,7 +414,8 @@ void Conv2d(
      ElementB, LayoutB,
      ElementC, LayoutC,
      ElementCompute,
-      ElementAccumulator, 
+      ElementAccumulator,
+      ElementD,
      ConvertOp, InnerProductOp
    >(problem_size, tensor_A, tensor_B, tensor_C, tensor_D, alpha, beta);
    break;
@@ -421,6 +427,7 @@ void Conv2d(
      ElementC, LayoutC,
      ElementCompute,
      ElementAccumulator,
+      ElementD,
      ConvertOp, InnerProductOp
    >(problem_size, tensor_A, tensor_B, tensor_C, tensor_D, alpha, beta);
    break;
@@ -431,7 +438,8 @@ void Conv2d(
      ElementB, LayoutB,
      ElementC, LayoutC,
      ElementCompute,
-      ElementAccumulator, 
+      ElementAccumulator,
+      ElementD,
      ConvertOp, InnerProductOp
    >(problem_size, tensor_A, tensor_B, tensor_C, tensor_D, alpha, beta);
    break;
--- a/tools/util/include/cutlass/util/reference/host/gemm_complex.h
+++ b/tools/util/include/cutlass/util/reference/host/gemm_complex.h
@@ -67,7 +67,8 @@ template <
  typename LayoutC,
  typename ScalarType,
  typename ComputeType,
-  typename ConvertOp = NumericConverter<ElementC, ScalarType>,
+  typename ElementD = ElementC,
+  typename ConvertOp = NumericConverter<ElementD, ScalarType>,
  typename InnerProductOp = multiply_add<ComputeType>
 >
 void GemmComplex(
@@ -79,7 +80,7 @@ void GemmComplex(
  ComplexTransform transform_b,
  ScalarType beta,
  TensorRef<ElementC, LayoutC> tensor_c,
-  TensorRef<ElementC, LayoutC> tensor_d,
+  TensorRef<ElementD, LayoutC> tensor_d,
  ComputeType initial_accum,
  int batch_count = 1,
  int64_t batch_stride_A = 0,
@@ -185,7 +186,8 @@ template <
  typename LayoutB,
  typename ElementC,
  typename LayoutC,
-  typename ScalarType
+  typename ScalarType,
+  typename ElementD = ElementC
 >
 void GemmComplex(
  gemm::GemmCoord problem_size,
@@ -196,7 +198,7 @@ void GemmComplex(
  ComplexTransform transform_b,
  ScalarType beta,
  TensorRef<ElementC, LayoutC> tensor_c,
-  TensorRef<ElementC, LayoutC> tensor_d) {
+  TensorRef<ElementD, LayoutC> tensor_d) {

  GemmComplex(problem_size, alpha, tensor_a, transform_a, tensor_b, transform_b, beta, tensor_c, tensor_d, ScalarType(0));
 }
--- a/tools/util/include/cutlass/util/reference/host/gett.hpp
+++ b/tools/util/include/cutlass/util/reference/host/gett.hpp
@@ -1,5 +1,5 @@
 /***************************************************************************************************
- * Copyright (c) 2017 - 2023 NVIDIA CORPORATION & AFFILIATES. All rights reserved.
+ * Copyright (c) 2023 - 2023 NVIDIA CORPORATION & AFFILIATES. All rights reserved.
 * SPDX-License-Identifier: BSD-3-Clause
 *
 * Redistribution and use in source and binary forms, with or without
@@ -38,6 +38,7 @@

 #include "cutlass/complex.h"
 #include "cutlass/numeric_conversion.h"
+#include "cutlass/epilogue/thread/activation.h"

 #include "cute/tensor.hpp"

@@ -75,7 +76,11 @@ template<
  class ElementAccumulator_,
  class ElementCompute_,
  class TensorC_,                                                                                          // (M, N, L)
-  class TensorD_                                                                                           // (M, N, L)
+  class TensorD_,                                                                                          // (M, N, L)
+  class TensorBias_,                                                                                       //    (M, 1)
+  class TensorT_,                                                                                          // (M, N, L)
+  class ActivationFunctor_ = cutlass::epilogue::thread::Identity<ElementCompute_>,
+  class BiasBinaryOp_ = cutlass::plus<ElementCompute_>
 >
 struct GettEpilogueParams {
  using ElementScalar = ElementScalar_;
@@ -83,15 +88,26 @@ struct GettEpilogueParams {
  using ElementCompute = ElementCompute_;
  using TensorC = TensorC_;
  using TensorD = TensorD_;
+  using TensorBias = TensorBias_;
+  using TensorT = TensorT_;
+  using ActivationFunctor = ActivationFunctor_;
+  using BiasBinaryOp = BiasBinaryOp_;
+
  using EngineC = typename TensorC::engine_type;
  using LayoutC = typename TensorC::layout_type;
  using EngineD =  typename TensorD::engine_type;
  using LayoutD = typename TensorD::layout_type;
+  using EngineBias =  typename TensorBias::engine_type;
+  using LayoutBias = typename TensorBias::layout_type;
+  using EngineT =  typename TensorT::engine_type;
+  using LayoutT = typename TensorT::layout_type;
  ElementScalar alpha = ElementScalar(1);
  ElementScalar beta = ElementScalar(0);

  TensorC C{};
  TensorD D{};
+  TensorBias Bias{};
+  TensorT T{};
 };

 /////////////////////////////////////////////////////////////////////////////////////////////////
@@ -204,19 +220,33 @@ void gett_epilogue(
  using ElementC = typename EpilogueParams::EngineC::value_type;

  using ElementD = typename EpilogueParams::EngineD::value_type;
+  using ElementBias = typename EpilogueParams::EngineBias::value_type;
+  using ElementT = typename EpilogueParams::EngineT::value_type;
+
  using ElementScalar = typename EpilogueParams::ElementScalar;
+  using ActivationFunctor = typename EpilogueParams::ActivationFunctor;
+  using BiasBinaryOp = typename EpilogueParams::BiasBinaryOp;
+
  // Input related converter
  NumericConverter<ElementCompute, ElementAccumulator> accumulator_converter;
  NumericConverter<ElementCompute, ElementC> source_converter;
+  NumericConverter<ElementCompute, ElementBias> bias_converter;

  // Scale related converter
  NumericConverter<ElementCompute, ElementScalar> scale_converter;
  // Output related converter
  NumericConverter<ElementD, ElementCompute> destination_converter;
+  NumericConverter<ElementT, ElementCompute> temporary_converter;
  // Epilogue operations
  multiply_add<ElementCompute, ElementCompute, ElementCompute> epilogue_fma;
  multiplies<ElementCompute> mul;

+  // Activation operation
+  ActivationFunctor activation;
+
+  // Bias binary operation
+  BiasBinaryOp bias_op;
+
  // Do conversion
  ElementCompute converted_alpha = scale_converter(epilogue_params.alpha);
  ElementCompute converted_beta = scale_converter(epilogue_params.beta);
@@ -225,10 +255,24 @@ void gett_epilogue(
      if (m + m_b < cute::size<0>(epilogue_params.D.layout()) && n + n_b < cute::size<1>(epilogue_params.D.layout())) {
        // Convert every type to ElementCompute first, do compute, convert to output type, write it out
        ElementCompute converted_acc = accumulator_converter(acc[m_b][n_b]);
-        ElementCompute converted_src = source_converter(epilogue_params.C(m + m_b, n + n_b, l));

-        ElementScalar output = epilogue_fma(converted_alpha, converted_acc, ElementCompute(0));
-        output = epilogue_fma(converted_beta, converted_src, output);
+        ElementCompute output = mul(converted_alpha, converted_acc);
+        if (epilogue_params.Bias.data()) {
+          ElementCompute converted_bias = bias_converter(epilogue_params.Bias(m + m_b));
+          output = bias_op(output, converted_bias);
+        }
+
+        if (epilogue_params.C.data()) {
+          ElementCompute converted_src = source_converter(epilogue_params.C(m + m_b, n + n_b, l));
+          output = epilogue_fma(converted_beta, converted_src, output);
+        }
+
+        if (epilogue_params.T.data()) {
+          // Store intermediate output
+          epilogue_params.T(m + m_b, n + n_b, l) = temporary_converter(output);
+        }
+
+        output = activation(output);

        epilogue_params.D(m + m_b, n + n_b, l) = destination_converter(output);
      }
@@ -238,6 +282,14 @@ void gett_epilogue(

 /////////////////////////////////////////////////////////////////////////////////////////////////

+template <class TensorType>
+auto make_layout_rank3(const TensorType& tensor) {
+  // append a batch mode of size 1 if we do not have tensors that are rank 3
+  return make_layout(
+      make_shape(get<0>(tensor.shape()), get<1>(tensor.shape()), Int<1>{}),
+      make_stride(get<0>(tensor.stride()), get<1>(tensor.stride()), int64_t(cosize(tensor.layout()))));
+}
+
 /// GEMM - General Matrix-Matrix contraction without conjugation options
 template <
  class MainloopParams,
@@ -254,26 +306,20 @@ void Gemm3x(
  static_assert(rank(typename MainloopParams::LayoutA{}) == rank(typename EpilogueParams::LayoutC{}));

  if constexpr (rank(typename MainloopParams::LayoutA{}) == 2) {
-    // append a batch mode of size 1 if we do not have tensors that are rank 3
-    Layout layout_A = make_layout(
-      make_shape(get<0>(mainloop_params.A.shape()), get<1>(mainloop_params.A.shape()), Int<1>{}),
-      make_stride(get<0>(mainloop_params.A.stride()), get<1>(mainloop_params.A.stride()), int64_t(cosize(mainloop_params.A.layout()))));
+    Layout layout_A = make_layout_rank3(mainloop_params.A);
+    Layout layout_B = make_layout_rank3(mainloop_params.B);
+    Layout layout_C = make_layout_rank3(epilogue_params.C);
+    Layout layout_D = make_layout_rank3(epilogue_params.D);
+    Layout layout_Bias = make_layout_rank3(epilogue_params.Bias);
+    Layout layout_T = make_layout_rank3(epilogue_params.T);

-    Layout layout_B = make_layout(
-      make_shape(get<0>(mainloop_params.B.shape()), get<1>(mainloop_params.B.shape()), Int<1>{}),
-      make_stride(get<0>(mainloop_params.B.stride()), get<1>(mainloop_params.B.stride()), int64_t(cosize(mainloop_params.B.layout()))));
-
-    Layout layout_C = make_layout(
-      make_shape(get<0>(epilogue_params.C.shape()), get<1>(epilogue_params.C.shape()), Int<1>{}),
-      make_stride(get<0>(epilogue_params.C.stride()), get<1>(epilogue_params.C.stride()), int64_t(cosize(epilogue_params.C.layout()))));
-
-    Layout layout_D = make_layout(
-      make_shape(get<0>(epilogue_params.D.shape()), get<1>(epilogue_params.D.shape()), Int<1>{}),
-      make_stride(get<0>(epilogue_params.D.stride()), get<1>(epilogue_params.D.stride()), int64_t(cosize(epilogue_params.D.layout()))));
    auto TensorA = make_tensor(mainloop_params.A.data(), layout_A);
    auto TensorB = make_tensor(mainloop_params.B.data(), layout_B);
    auto TensorC = make_tensor(epilogue_params.C.data(), layout_C);
    auto TensorD = make_tensor(epilogue_params.D.data(), layout_D);
+    auto TensorBias = make_tensor(epilogue_params.Bias.data(), layout_Bias);
+    auto TensorT = make_tensor(epilogue_params.T.data(), layout_T);
+
    // Reconstruct mainloop params
    GettMainloopParams<typename MainloopParams::ElementAccumulator,
                       decltype(TensorA),
@@ -288,12 +334,16 @@ void Gemm3x(
                       typename EpilogueParams::ElementAccumulator,
                       typename EpilogueParams::ElementCompute,
                       decltype(TensorC),
-                       decltype(TensorD)
+                       decltype(TensorD),
+                       decltype(TensorBias),
+                       decltype(TensorT)
                      >
        epilogue_params_converted{epilogue_params.alpha,
                                  epilogue_params.beta,
                                  TensorC,
-                                  TensorD
+                                  TensorD,
+                                  TensorBias,
+                                  TensorT
                                  };

    Gett(mainloop_params_converted, epilogue_params_converted);
--- a/tools/util/include/cutlass/util/reference/host/tensor_compare.h
+++ b/tools/util/include/cutlass/util/reference/host/tensor_compare.h
@@ -39,11 +39,11 @@

 // Cutlass includes
 #include "cutlass/cutlass.h"
+#include "cutlass/relatively_equal.h"
 #include "cutlass/tensor_view.h"
 #include "cutlass/tensor_view_planar_complex.h"

 #include "cutlass/util/distribution.h"
-//#include "cutlass/util/type_traits.h"
 #include "tensor_foreach.h"

 namespace cutlass {
@@ -83,10 +83,55 @@ struct TensorEqualsFunc {

    Element lhs_ = lhs.at(coord);
    Element rhs_ = rhs.at(coord);
-    
+
    if (lhs_ != rhs_) {
      result = false;
-    }  
+    }
+  }
+
+  /// Returns true if equal
+  operator bool() const {
+    return result;
+  }
+};
+
+template <
+  typename Element,               ///< Element type
+  typename Layout>                ///< Layout function
+struct TensorRelativelyEqualsFunc {
+
+  //
+  // Data members
+  //
+
+  TensorView<Element, Layout> lhs;
+  TensorView<Element, Layout> rhs;
+  Element epsilon;
+  Element nonzero_floor;
+  bool result;
+
+  /// Ctor
+  TensorRelativelyEqualsFunc(
+    TensorView<Element, Layout> const &lhs_,
+    TensorView<Element, Layout> const &rhs_,
+    Element epsilon_,
+    Element nonzero_floor_
+  ) :
+    lhs(lhs_),
+    rhs(rhs_),
+    epsilon(epsilon_),
+    nonzero_floor(nonzero_floor_),
+    result(true) { }
+
+  /// Visits a coordinate
+  void operator()(Coord<Layout::kRank> const &coord) {
+
+    Element lhs_ = lhs.at(coord);
+    Element rhs_ = rhs.at(coord);
+
+    if (!relatively_equal(lhs_, rhs_, epsilon, nonzero_floor)) {
+      result = false;
+    }
  }

  /// Returns true if equal
@@ -104,7 +149,7 @@ template <
  typename Element,               ///< Element type
  typename Layout>                ///< Layout function
 bool TensorEquals(
-  TensorView<Element, Layout> const &lhs, 
+  TensorView<Element, Layout> const &lhs,
  TensorView<Element, Layout> const &rhs) {

  // Extents must be identical
@@ -126,7 +171,7 @@ template <
  typename Element,               ///< Element type
  typename Layout>                ///< Layout function
 bool TensorEquals(
-  TensorViewPlanarComplex<Element, Layout> const &lhs, 
+  TensorViewPlanarComplex<Element, Layout> const &lhs,
  TensorViewPlanarComplex<Element, Layout> const &rhs) {

  // Extents must be identical
@@ -135,7 +180,7 @@ bool TensorEquals(
  }

  detail::TensorEqualsFunc<Element, Layout> real_func(
-    {lhs.data(), lhs.layout(), lhs.extent()}, 
+    {lhs.data(), lhs.layout(), lhs.extent()},
    {rhs.data(), rhs.layout(), rhs.extent()}
  );

@@ -164,12 +209,85 @@ bool TensorEquals(
 ///////////////////////////////////////////////////////////////////////////////////////////////////
 ///////////////////////////////////////////////////////////////////////////////////////////////////

+/// Returns true if two tensor views are relatively equal.
+template <
+  typename Element,               ///< Element type
+  typename Layout>                ///< Layout function
+bool TensorRelativelyEquals(
+  TensorView<Element, Layout> const &lhs,
+  TensorView<Element, Layout> const &rhs,
+  Element epsilon,
+  Element nonzero_floor) {
+
+  // Extents must be identical
+  if (lhs.extent() != rhs.extent()) {
+    return false;
+  }
+
+  detail::TensorRelativelyEqualsFunc<Element, Layout> func(lhs, rhs, epsilon, nonzero_floor);
+  TensorForEach(
+    lhs.extent(),
+    func
+  );
+
+  return bool(func);
+}
+
+/// Returns true if two tensor views are relatively equal.
+template <
+  typename Element,               ///< Element type
+  typename Layout>                ///< Layout function
+bool TensorRelativelyEquals(
+  TensorViewPlanarComplex<Element, Layout> const &lhs,
+  TensorViewPlanarComplex<Element, Layout> const &rhs,
+  Element epsilon,
+  Element nonzero_floor) {
+
+  // Extents must be identical
+  if (lhs.extent() != rhs.extent()) {
+    return false;
+  }
+
+  detail::TensorRelativelyEqualsFunc<Element, Layout> real_func(
+    {lhs.data(), lhs.layout(), lhs.extent()},
+    {rhs.data(), rhs.layout(), rhs.extent()},
+    epsilon,
+    nonzero_floor
+  );
+
+  TensorForEach(
+    lhs.extent(),
+    real_func
+  );
+
+  if (!bool(real_func)) {
+    return false;
+  }
+
+  detail::TensorEqualsFunc<Element, Layout> imag_func(
+    {lhs.data() + lhs.imaginary_stride(), lhs.layout(), lhs.extent()},
+    {rhs.data() + rhs.imaginary_stride(), rhs.layout(), rhs.extent()},
+    epsilon,
+    nonzero_floor
+  );
+
+  TensorForEach(
+    lhs.extent(),
+    imag_func
+  );
+
+  return bool(imag_func);
+}
+
+///////////////////////////////////////////////////////////////////////////////////////////////////
+///////////////////////////////////////////////////////////////////////////////////////////////////
+
 /// Returns true if two tensor views are NOT equal.
 template <
  typename Element,               ///< Element type
  typename Layout>                ///< Layout function
 bool TensorNotEquals(
-  TensorView<Element, Layout> const &lhs, 
+  TensorView<Element, Layout> const &lhs,
  TensorView<Element, Layout> const &rhs) {

  // Extents must be identical
@@ -191,7 +309,7 @@ template <
  typename Element,               ///< Element type
  typename Layout>                ///< Layout function
 bool TensorNotEquals(
-  TensorViewPlanarComplex<Element, Layout> const &lhs, 
+  TensorViewPlanarComplex<Element, Layout> const &lhs,
  TensorViewPlanarComplex<Element, Layout> const &rhs) {

  return !TensorEquals(lhs, rhs);
@@ -235,7 +353,7 @@ struct TensorContainsFunc {

    if (view.at(coord) == value) {
      if (!contains) {
-        location = coord; 
+        location = coord;
      }
      contains = true;
    }
--- a/tools/util/include/cutlass/util/tensor_view_io.h
+++ b/tools/util/include/cutlass/util/tensor_view_io.h
@@ -96,12 +96,16 @@ inline std::ostream & TensorView_WriteRank(

    if (rank + 2 == Layout::kRank) {
      // Write least significant ranks asa matrix with rows delimited by "\n"
-      out << (idx ? ",\n" : "");
+      if (idx) {
+        out << ",\n";
+      }
      TensorView_WriteLeastSignificantRank(out, view, coord, rank + 1, width);
    }
    else {
      // Higher ranks are separated by newlines
-      out << (idx ? ",\n\n" : "");
+      if (idx) {
+        out << ",\n\n";
+      }
      TensorView_WriteRank(out, view, coord, rank + 1, width);
    }
  }
@@ -166,12 +170,16 @@ inline std::ostream & TensorViewPlanarComplex_WriteRank(

    if (rank + 2 == Layout::kRank) {
      // Write least significant ranks asa matrix with rows delimited by ";\n"
-      out << (idx ? ";\n" : "");
+      if (idx) {
+        out << ";\n";
+      }
      TensorViewPlanarComplex_WriteLeastSignificantRank(out, view, coord, rank + 1, width);
    }
    else {
      // Higher ranks are separated by newlines
-      out << (idx ? "\n" : "");
+      if (idx) {
+        out << "\n";
+      }
      TensorViewPlanarComplex_WriteRank(out, view, coord, rank + 1, width);
    }
  }