#include <vtkm/cont/ArrayHandleCast.h>
#include <vtkm/cont/ArrayHandleGroupVecVariable.h>
#include <vtkm/cont/CellSetExplicit.h>

#include <vtkm/exec/CellFace.h>

#include <vtkm/Hash.h>

#include <vtkm/worklet/AverageByKey.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
#include <vtkm/worklet/DispatcherReduceByKey.h>
#include <vtkm/worklet/Keys.h>
#include <vtkm/worklet/ScatterCounting.h>

#include <vtkm/filter/FilterDataSet.h>

#include <vtkm/cont/testing/MakeTestDataSet.h>
#include <vtkm/cont/testing/Testing.h>

//#define CHECK_COLLISIONS

namespace vtkm {
namespace worklet {

namespace {


struct ExtractFaces
{
  struct CountFaces : vtkm::worklet::WorkletMapPointToCell
  {
    typedef void ControlSignature(CellSetIn cellSet, FieldOut<> numFaces);
    typedef _2 ExecutionSignature(CellShape);
    using InputDomain = _1;

    template<typename CellShapeTag>
    VTKM_EXEC_CONT
    vtkm::IdComponent operator()(CellShapeTag cellShape) const
    {
      return vtkm::exec::CellFaceNumberOfFaces(cellShape, *this);
    }
  };

  class FaceHashes : public vtkm::worklet::WorkletMapPointToCell
  {
  public:
    typedef void ControlSignature(CellSetIn cellSet,
                                  FieldOut<> hashValues);
    typedef _2 ExecutionSignature(CellShape cellShape,
                                  PointIndices globalPointIndices,
                                  VisitIndex localFaceIndex);
    using InputDomain = _1;

    using ScatterType = vtkm::worklet::ScatterCounting;

    template <typename CellShapeTag, typename PointIndexVecType>
    VTKM_EXEC
    vtkm::HashType operator()(CellShapeTag cellShape,
                              const PointIndexVecType& globalPointIndicesForCell,
                              vtkm::IdComponent localFaceIndex) const
    {
      //// PAUSE-EXAMPLE
#ifndef CHECK_COLLISIONS
      //// RESUME-EXAMPLE
      return vtkm::Hash(
            vtkm::exec::CellFaceCanonicalId(localFaceIndex,
                                            cellShape,
                                            globalPointIndicesForCell,
                                            *this));
      //// PAUSE-EXAMPLE
#else // ! CHECK_COLLISIONS
      // Intentionally use a bad hash value to cause collisions to check to make
      // sure that collision resolution works.
      return vtkm::HashType(
            vtkm::exec::CellFaceCanonicalId(localFaceIndex,
                                            cellShape,
                                            globalPointIndicesForCell,
                                            *this)[0]);
#endif // !CHECK_COLLISIONS
      //// RESUME-EXAMPLE
    }
  };

  class FaceHashCollisions : public vtkm::worklet::WorkletReduceByKey
  {
  public:
    typedef void ControlSignature(KeysIn keys,
                                  WholeCellSetIn<> inputCells,
                                  ValuesIn<> originCells,
                                  ValuesIn<> originFaces,
                                  ValuesOut<> localFaceIndices,
                                  ReducedValuesOut<> numFaces);
    typedef _6 ExecutionSignature(_2 inputCells,
                                  _3 originCells,
                                  _4 originFaces,
                                  _5 localFaceIndices);
    using InputDomain = _1;

    template <typename CellSetType,
              typename OriginCellsType,
              typename OriginFacesType,
              typename localFaceIndicesType>
    VTKM_EXEC
    vtkm::IdComponent operator()(const CellSetType& cellSet,
                                 const OriginCellsType& originCells,
                                 const OriginFacesType& originFaces,
                                 localFaceIndicesType& localFaceIndices) const
    {
      vtkm::IdComponent numFacesInHash = localFaceIndices.GetNumberOfComponents();

      // Sanity checks.
      VTKM_ASSERT(originCells.GetNumberOfComponents() == numFacesInHash);
      VTKM_ASSERT(originFaces.GetNumberOfComponents() == numFacesInHash);

      // Clear out localFaceIndices
      for (vtkm::IdComponent index = 0; index < numFacesInHash; ++index)
      {
        localFaceIndices[index] = -1;
      }

      // Count how many unique faces there are and create an id for each;
      vtkm::IdComponent numUniqueFaces = 0;
      for (vtkm::IdComponent firstFaceIndex = 0;
           firstFaceIndex < numFacesInHash;
           ++firstFaceIndex)
      {
        if (localFaceIndices[firstFaceIndex] == -1)
        {
          vtkm::IdComponent faceId = numUniqueFaces;
          localFaceIndices[firstFaceIndex] = faceId;
          // Find all matching faces.
          vtkm::Id firstCellIndex = originCells[firstFaceIndex];
          vtkm::Id3 canonicalFaceId =
              vtkm::exec::CellFaceCanonicalId(
                originFaces[firstFaceIndex],
                cellSet.GetCellShape(firstCellIndex),
                cellSet.GetIndices(firstCellIndex),
                *this);
          for (vtkm::IdComponent laterFaceIndex = firstFaceIndex+1;
               laterFaceIndex < numFacesInHash;
               ++laterFaceIndex)
          {
            vtkm::Id laterCellIndex = originCells[laterFaceIndex];
            vtkm::Id3 otherCanonicalFaceId =
                vtkm::exec::CellFaceCanonicalId(
                  originFaces[laterFaceIndex],
                  cellSet.GetCellShape(laterCellIndex),
                  cellSet.GetIndices(laterCellIndex),
                  *this);
            if (canonicalFaceId == otherCanonicalFaceId)
            {
              localFaceIndices[laterFaceIndex] = faceId;
            }
          }
          ++numUniqueFaces;
        }
      }

      return numUniqueFaces;
    }
  };

  ////
  //// BEGIN-EXAMPLE GenerateMeshVariableShapeCountPointsInFace.cxx
  ////
  class CountPointsInFace : public vtkm::worklet::WorkletReduceByKey
  {
  public:
    typedef void ControlSignature(KeysIn keys,
                                  WholeCellSetIn<> inputCells,
                                  ValuesIn<> originCells,
                                  ValuesIn<> originFaces,
                                  ValuesIn<> localFaceIndices,
                                  ReducedValuesOut<> faceShape,
                                  ReducedValuesOut<> numPointsInEachFace);
    typedef void ExecutionSignature(_2 inputCells,
                                    _3 originCell,
                                    _4 originFace,
                                    _5 localFaceIndices,
                                    VisitIndex localFaceIndex,
                                    _6 faceShape,
                                    _7 numPointsInFace);
    using InputDomain = _1;

    using ScatterType = vtkm::worklet::ScatterCounting;

    template <typename CellSetType,
              typename OriginCellsType,
              typename OriginFacesType,
              typename LocalFaceIndicesType>
    VTKM_EXEC
    void operator()(const CellSetType& cellSet,
                    const OriginCellsType& originCells,
                    const OriginFacesType& originFaces,
                    const LocalFaceIndicesType& localFaceIndices,
                    vtkm::IdComponent localFaceIndex,
                    vtkm::UInt8& faceShape,
                    vtkm::IdComponent& numPointsInFace) const
    {
      // Find the first face that matches the index given.
      for (vtkm::IdComponent faceIndex = 0; ; ++faceIndex)
      {
        if (localFaceIndices[faceIndex] == localFaceIndex)
        {
          vtkm::Id cellIndex = originCells[faceIndex];
          faceShape =
              vtkm::exec::CellFaceShape(originFaces[faceIndex],
                                        cellSet.GetCellShape(cellIndex),
                                        *this);
          numPointsInFace =
              vtkm::exec::CellFaceNumberOfPoints(originFaces[faceIndex],
                                                 cellSet.GetCellShape(cellIndex),
                                                 *this);
          break;
        }
      }
    }
  };
  ////
  //// END-EXAMPLE GenerateMeshVariableShapeCountPointsInFace.cxx
  ////

  ////
  //// BEGIN-EXAMPLE GenerateMeshVariableShapeGenIndices.cxx
  ////
  class FaceIndices : public vtkm::worklet::WorkletReduceByKey
  {
  public:
    typedef void ControlSignature(KeysIn keys,
                                  WholeCellSetIn<> inputCells,
                                  ValuesIn<> originCells,
                                  ValuesIn<> originFaces,
                                  ValuesIn<> localFaceIndices,
                                  ReducedValuesOut<> connectivityOut);
    typedef void ExecutionSignature(_2 inputCells,
                                    _3 originCell,
                                    _4 originFace,
                                    _5 localFaceIndices,
                                    VisitIndex localFaceIndex,
                                    _6 connectivityOut);
    using InputDomain = _1;

    using ScatterType = vtkm::worklet::ScatterCounting;

    template <typename CellSetType,
              typename OriginCellsType,
              typename OriginFacesType,
              typename LocalFaceIndicesType,
              typename ConnectivityVecType>
    VTKM_EXEC
    void operator()(const CellSetType& cellSet,
                    const OriginCellsType& originCells,
                    const OriginFacesType& originFaces,
                    const LocalFaceIndicesType& localFaceIndices,
                    vtkm::IdComponent localFaceIndex,
                    ConnectivityVecType& connectivityOut) const
    {
      // Find the first face that matches the index given and return it.
      for (vtkm::IdComponent faceIndex = 0; ; ++faceIndex)
      {
        if (localFaceIndices[faceIndex] == localFaceIndex)
        {
          vtkm::Id cellIndex = originCells[faceIndex];
          vtkm::IdComponent faceInCellIndex = originFaces[faceIndex];
          auto cellShape = cellSet.GetCellShape(cellIndex);
          vtkm::IdComponent numPointsInFace =
              connectivityOut.GetNumberOfComponents();

          VTKM_ASSERT(numPointsInFace == vtkm::exec::CellFaceNumberOfPoints(
                        faceInCellIndex, cellShape, *this));

          auto globalPointIndicesForCell = cellSet.GetIndices(cellIndex);
          for (vtkm::IdComponent localPointI = 0;
               localPointI < numPointsInFace;
               ++localPointI)
          {
            vtkm::IdComponent pointInCellIndex =
                vtkm::exec::CellFaceLocalIndex(localPointI,
                                               faceInCellIndex,
                                               cellShape,
                                               *this);
            connectivityOut[localPointI] =
                globalPointIndicesForCell[pointInCellIndex];
          }

          break;
        }
      }
    }
  };
  ////
  //// END-EXAMPLE GenerateMeshVariableShapeGenIndices.cxx
  ////

  ////
  //// BEGIN-EXAMPLE GenerateMeshVariableShapeInvoke.cxx
  ////
  template<typename InCellSetType, typename Device>
  VTKM_CONT
  vtkm::cont::CellSetExplicit<> Run(const InCellSetType& inCellSet, Device)
  {
    VTKM_IS_DYNAMIC_OR_STATIC_CELL_SET(InCellSetType);

    using OutCellSetType = vtkm::cont::CellSetExplicit<>;

    // First, count the faces in each cell.
    vtkm::cont::ArrayHandle<vtkm::IdComponent> faceCounts;
    vtkm::worklet::DispatcherMapTopology<CountFaces, Device> countFaceDispatcher;
    countFaceDispatcher.Invoke(inCellSet, faceCounts);

    vtkm::worklet::ScatterCounting scatter(faceCounts, Device());
    this->OutputToInputCellMap =
        scatter.GetOutputToInputMap(inCellSet.GetNumberOfCells());
    vtkm::worklet::ScatterCounting::VisitArrayType outputToInputFaceMap =
        scatter.GetVisitArray(inCellSet.GetNumberOfCells());

    // Second, for each face, extract a hash.
    vtkm::cont::ArrayHandle<vtkm::HashType> hashValues;

    vtkm::worklet::DispatcherMapTopology<FaceHashes, Device>
        FaceHashesDispatcher(scatter);
    FaceHashesDispatcher.Invoke(inCellSet, hashValues);

    // Third, use a Keys object to combine all like hashes.
    this->CellToFaceKeys = vtkm::worklet::Keys<vtkm::HashType>(hashValues, Device());

    // Fourth, use a reduce-by-key to collect like hash values, resolve collisions,
    // and count the number of unique faces associated with each hash.
    vtkm::cont::ArrayHandle<vtkm::IdComponent> numUniqueFacesInEachHash;

    vtkm::worklet::DispatcherReduceByKey<FaceHashCollisions>
        faceHashCollisionDispatcher;
    faceHashCollisionDispatcher.Invoke(this->CellToFaceKeys,
                                       inCellSet,
                                       this->OutputToInputCellMap,
                                       outputToInputFaceMap,
                                       this->LocalFaceIndices,
                                       numUniqueFacesInEachHash);

    // Fifth, use a reduce-by-key to count the number of points in each unique face.
    // Also identify the shape of each face.
    this->HashCollisionScatter.reset(
          new vtkm::worklet::ScatterCounting(numUniqueFacesInEachHash, Device()));

    OutCellSetType::ShapeArrayType shapeArray;
    OutCellSetType::NumIndicesArrayType numPointsInEachFace;

    vtkm::worklet::DispatcherReduceByKey<CountPointsInFace, Device>
        countPointsDispatcher(*this->HashCollisionScatter);
    countPointsDispatcher.Invoke(this->CellToFaceKeys,
                                 inCellSet,
                                 this->OutputToInputCellMap,
                                 outputToInputFaceMap,
                                 this->LocalFaceIndices,
                                 shapeArray,
                                 numPointsInEachFace);

    // Sixth, convert the numPointsInEachFace array to an offsets array and use that
    // to create an ArrayHandleGroupVecVariable.
    ////
    //// BEGIN-EXAMPLE GenerateMeshVariableShapeOffsetsArray.cxx
    ////
    vtkm::cont::ArrayHandle<vtkm::Id> offsetsArray;
    vtkm::Id connectivityArraySize;
    vtkm::cont::ConvertNumComponentsToOffsets(
          numPointsInEachFace, offsetsArray, connectivityArraySize, Device());
    OutCellSetType::ConnectivityArrayType connectivityArray;
    connectivityArray.Allocate(connectivityArraySize);
    auto connectivityArrayVecs = vtkm::cont::make_ArrayHandleGroupVecVariable(
          connectivityArray, offsetsArray);
    ////
    //// END-EXAMPLE GenerateMeshVariableShapeOffsetsArray.cxx
    ////

    // Seventh, use a reduce-by-key to extract indices for each unique face.
    vtkm::worklet::DispatcherReduceByKey<FaceIndices, Device>
        faceIndicesDispatcher(*this->HashCollisionScatter);
    faceIndicesDispatcher.Invoke(this->CellToFaceKeys,
                                 inCellSet,
                                 this->OutputToInputCellMap,
                                 outputToInputFaceMap,
                                 this->LocalFaceIndices,
                                 connectivityArrayVecs);

    // Eigth, use the created connectivity array and others to build a cell set.
    OutCellSetType outCellSet(inCellSet.GetName());
    outCellSet.Fill(inCellSet.GetNumberOfPoints(),
                    shapeArray,
                    numPointsInEachFace,
                    connectivityArray);

    return outCellSet;
  }
  ////
  //// END-EXAMPLE GenerateMeshVariableShapeInvoke.cxx
  ////

  class AverageCellField : public vtkm::worklet::WorkletReduceByKey
  {
  public:
    typedef void ControlSignature(KeysIn keys,
                                  ValuesIn<> inFieldValues,
                                  ValuesIn<> localFaceIndices,
                                  ReducedValuesOut<> averagedField);
    typedef _4 ExecutionSignature(_2 inFieldValues,
                                  _3 localFaceIndices,
                                  VisitIndex localFaceIndex);
    using InputDomain = _1;

    using ScatterType = vtkm::worklet::ScatterCounting;

    template <typename InFieldValuesType, typename LocalFaceIndicesType>
    VTKM_EXEC
    typename InFieldValuesType::ComponentType
    operator()(const InFieldValuesType& inFieldValues,
               const LocalFaceIndicesType& localFaceIndices,
               vtkm::IdComponent localFaceIndex) const
    {
      using FieldType = typename InFieldValuesType::ComponentType;

      FieldType averageField = FieldType(0);
      vtkm::IdComponent numValues = 0;
      for (vtkm::IdComponent reduceIndex = 0;
           reduceIndex < inFieldValues.GetNumberOfComponents();
           ++reduceIndex)
      {
        if (localFaceIndices[reduceIndex] == localFaceIndex)
        {
          FieldType fieldValue = inFieldValues[reduceIndex];
          averageField = averageField + fieldValue;
          ++numValues;
        }
      }
      VTKM_ASSERT(numValues > 0);
      return averageField / numValues;
    }
  };

  ////
  //// BEGIN-EXAMPLE GenerateMeshVariableShapeMapCellField.cxx
  ////
  template <typename ValueType, typename Storage, typename Device>
  VTKM_CONT
  vtkm::cont::ArrayHandle<ValueType> ProcessCellField(
      const vtkm::cont::ArrayHandle<ValueType, Storage>& inCellField,
      Device) const
  {
    vtkm::cont::ArrayHandle<ValueType> averageField;
    vtkm::worklet::DispatcherReduceByKey<AverageCellField, Device>
        dispatcher(*this->HashCollisionScatter);
    dispatcher.Invoke(
          this->CellToFaceKeys,
          vtkm::cont::make_ArrayHandlePermutation(this->OutputToInputCellMap,
                                                  inCellField),
          this->LocalFaceIndices,
          averageField);
    return averageField;
  }
  ////
  //// END-EXAMPLE GenerateMeshVariableShapeMapCellField.cxx
  ////

  vtkm::worklet::ScatterCounting::OutputToInputMapType OutputToInputCellMap;
  vtkm::worklet::Keys<vtkm::HashType> CellToFaceKeys;
  vtkm::cont::ArrayHandle<vtkm::IdComponent> LocalFaceIndices;
  std::shared_ptr<vtkm::worklet::ScatterCounting> HashCollisionScatter;
};

} // anonymous namespace

}
} // namespace vtkm::worklet

namespace vtkm {
namespace filter {

//// PAUSE-EXAMPLE
namespace {

//// RESUME-EXAMPLE
class ExtractFaces : public vtkm::filter::FilterDataSet<ExtractFaces>
{
public:
  template <typename Policy, typename Device>
  VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& inData,
                                          vtkm::filter::PolicyBase<Policy> policy,
                                          Device);

  template <typename T, typename StorageType, typename Policy, typename Device>
  VTKM_CONT bool DoMapField(vtkm::cont::DataSet& result,
                            const vtkm::cont::ArrayHandle<T, StorageType>& input,
                            const vtkm::filter::FieldMetadata& fieldMeta,
                            const vtkm::filter::PolicyBase<Policy>& policy,
                            Device);

private:
  vtkm::worklet::ExtractFaces Worklet;
};

//// PAUSE-EXAMPLE
} // anonymous namespace
//// RESUME-EXAMPLE
}
} // namespace vtkm::filter


namespace vtkm {
namespace filter {

//// PAUSE-EXAMPLE
namespace {

//// RESUME-EXAMPLE
template <typename Policy, typename Device>
inline VTKM_CONT vtkm::cont::DataSet
ExtractFaces::DoExecute(const vtkm::cont::DataSet& inData,
                        vtkm::filter::PolicyBase<Policy> policy,
                        Device)
{
  VTKM_IS_DEVICE_ADAPTER_TAG(Device);

  const vtkm::cont::DynamicCellSet& inCells =
      inData.GetCellSet(this->GetActiveCellSetIndex());

  vtkm::cont::CellSetExplicit<> outCells =
      this->Worklet.Run(vtkm::filter::ApplyPolicy(inCells, policy), Device());

  vtkm::cont::DataSet outData;

  outData.AddCellSet(outCells);

  for (vtkm::IdComponent coordSystemIndex = 0;
       coordSystemIndex < inData.GetNumberOfCoordinateSystems();
       ++coordSystemIndex)
  {
    outData.AddCoordinateSystem(inData.GetCoordinateSystem(coordSystemIndex));
  }

  return outData;
}

template <typename T, typename StorageType, typename Policy, typename Device>
inline VTKM_CONT
bool ExtractFaces::DoMapField(vtkm::cont::DataSet& result,
                              const vtkm::cont::ArrayHandle<T, StorageType>& input,
                              const vtkm::filter::FieldMetadata& fieldMeta,
                              const vtkm::filter::PolicyBase<Policy>&,
                              Device)
{
  vtkm::cont::Field output;

  if (fieldMeta.IsPointField())
  {
    output = fieldMeta.AsField(input); // pass through
  }
  else if (fieldMeta.IsCellField())
  {
    output = fieldMeta.AsField(this->Worklet.ProcessCellField(input, Device()));
  }
  else
  {
    return false;
  }

  result.AddField(output);

  return true;
}

//// PAUSE-EXAMPLE
} // anonymous namespace

//// RESUME-EXAMPLE
}
} // namespace vtkm::filter

namespace {

template <typename FaceVecLikeType>
bool IsFace(const FaceVecLikeType& faceIndices,
            const vtkm::Id3& expectedIndices)
{
  for (vtkm::IdComponent expectedIndex = 0; expectedIndex < 3; ++expectedIndex)
  {
    bool foundIndex = false;
    for (vtkm::IdComponent faceIndex = 0;
         faceIndex < faceIndices.GetNumberOfComponents();
         ++faceIndex)
    {
      if (expectedIndices[expectedIndex] == faceIndices[faceIndex])
      {
        foundIndex = true;
        break;
      }
    }
    if (!foundIndex) { return false; }
  }
  return true;
}

template <typename NumPointsInFaceType,
          typename ConnectivityPortalType,
          typename OffsetsPortalType>
vtkm::Id FindFace(const NumPointsInFaceType& numPointsInFace,
                  const ConnectivityPortalType& connectivity,
                  const OffsetsPortalType& offsets,
                  const vtkm::Id3& face)
{
  vtkm::Id faceIndex = -1;
  bool foundFace = false;

  for (vtkm::Id offsetIndex = 0;
       offsetIndex < offsets.GetNumberOfValues();
       ++offsetIndex)
  {
    if (IsFace(vtkm::VecFromPortal<ConnectivityPortalType>(
                 connectivity,
                 numPointsInFace.Get(offsetIndex),
                 offsets.Get(offsetIndex)),
               face))
    {
      faceIndex = offsetIndex;
      break;
    }
  }

  VTKM_TEST_ASSERT(faceIndex >= 0, "Did not find expected face.");

  for (vtkm::Id offsetIndex = faceIndex+1;
       offsetIndex < offsets.GetNumberOfValues();
       ++offsetIndex)
  {
    if (IsFace(vtkm::VecFromPortal<ConnectivityPortalType>(
                 connectivity,
                 numPointsInFace.Get(offsetIndex),
                 offsets.Get(offsetIndex)),
               face))
    {
      VTKM_TEST_FAIL("Face duplicated.");
    }
  }

  return faceIndex;
}

template <typename NumPointsInFaceType,
          typename ConnectivityPortalType,
          typename OffsetsPortalType,
          typename FieldPortalType>
void CheckFace(const NumPointsInFaceType& numPointsInFace,
               const ConnectivityPortalType& connectivity,
               const OffsetsPortalType& offsets,
               const FieldPortalType& field,
               const vtkm::Id3& face,
               const vtkm::Float32 expectedFieldValue)
{
  std::cout << "  Checking for face " << face << " with field value "
            << expectedFieldValue << std::endl;

  vtkm::Id faceIndex = FindFace(numPointsInFace, connectivity, offsets, face);

  vtkm::Float32 fieldValue = field.Get(faceIndex);
  VTKM_TEST_ASSERT(test_equal(expectedFieldValue, fieldValue), "Bad field value.");
}

void CheckOutput(const vtkm::cont::CellSetExplicit<>& cellSet,
                 const vtkm::cont::ArrayHandle<vtkm::Float32>& cellField)
{
  std::cout << "Num cells: " << cellSet.GetNumberOfCells() << std::endl;
  VTKM_TEST_ASSERT(cellSet.GetNumberOfCells() == 16, "Wrong # of cells.");

  auto numPointsPerFace = cellSet.GetNumIndicesArray(vtkm::TopologyElementTagPoint(),
                                                     vtkm::TopologyElementTagCell());
  auto connectivity = cellSet.GetConnectivityArray(vtkm::TopologyElementTagPoint(),
                                                   vtkm::TopologyElementTagCell());
  auto offsets = cellSet.GetIndexOffsetArray(vtkm::TopologyElementTagPoint(),
                                             vtkm::TopologyElementTagCell());

  std::cout << "Offsets:" << std::endl;
  vtkm::cont::printSummary_ArrayHandle(offsets, std::cout, true);
  std::cout << "Connectivity:" << std::endl;
  vtkm::cont::printSummary_ArrayHandle(connectivity, std::cout, true);

  std::cout << "Cell field:" << std::endl;
  vtkm::cont::printSummary_ArrayHandle(cellField, std::cout, true);

  auto nPointsP = numPointsPerFace.GetPortalConstControl();
  auto connectivityP = connectivity.GetPortalConstControl();
  auto offsetsP = offsets.GetPortalConstControl();
  auto fieldP = cellField.GetPortalConstControl();
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(0, 1, 2), 100.1f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(0, 1, 5), 100.1f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(0, 3, 7), 100.1f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(6, 2, 3), 100.1f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(6, 2, 1), 105.05f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(6, 5, 4), 115.3f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(8, 1, 2), 110.0f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(8, 1, 5), 110.0f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(8, 6, 2), 110.0f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(8, 6, 5), 115.1f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(10, 6, 8), 120.2f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(10, 5, 8), 120.2f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(10, 6, 5), 125.35f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(9, 4, 7), 130.5f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(9, 4, 5), 130.5f);
  CheckFace(nPointsP, connectivityP, offsetsP, fieldP, vtkm::Id3(9, 4, 10), 130.5f);
}

//void TryWorklet()
//{
//  std::cout << std::endl << "Trying calling worklet." << std::endl;
//  vtkm::cont::DataSet inDataSet =
//      vtkm::cont::testing::MakeTestDataSet().Make3DExplicitDataSet5();
//  vtkm::cont::CellSetExplicit<> inCellSet;
//  inDataSet.GetCellSet().CopyTo(inCellSet);

//  vtkm::worklet::ExtractFaces worklet;
//  vtkm::cont::CellSetExplicit<> outCellSet =
//      worklet.Run(inCellSet, VTKM_DEFAULT_DEVICE_ADAPTER_TAG());
//  CheckOutput(outCellSet);
//}

void TryFilter()
{
  std::cout << std::endl << "Trying calling filter." << std::endl;
  vtkm::cont::DataSet inDataSet =
      vtkm::cont::testing::MakeTestDataSet().Make3DExplicitDataSet5();

  vtkm::filter::ExtractFaces filter;

  // NOTE 2018-03-21: I expect this to fail in the short term. Right now no fields
  // are copied from input to output. The default should be changed to copy them
  // all.
  vtkm::cont::DataSet outDataSet = filter.Execute(inDataSet);

  vtkm::cont::CellSetExplicit<> outCellSet;
  outDataSet.GetCellSet().CopyTo(outCellSet);

  vtkm::cont::Field outCellField = outDataSet.GetField("cellvar");
  vtkm::cont::ArrayHandle<vtkm::Float32> outCellData;
  outCellField.GetData().CopyTo(outCellData);

  CheckOutput(outCellSet, outCellData);
}

void DoTest()
{
//  TryWorklet();
  TryFilter();
}

} // anonymous namespace

int GenerateMeshVariableShape(int, char*[])
{
  return vtkm::cont::testing::Testing::Run(DoTest);
}
