Merge topic 'zfp'

865942f18 Home stretch of warnings...(?)
f5569fb93 More warnings.
6beefb4ab Technically need default initialization.
6f86d21e3 clang is my nemesis.
ce00cc6e6 straggling function variables.
3dae0b08f CUDA can still be finicky.
e7c1151c7 So many warnings.
f5ca28f5c More compiler warnings.
...

Acked-by: Kitware Robot <kwrobot@kitware.com>
Acked-by: Kenneth Moreland <kmorel@sandia.gov>
Merge-request: !1482
This commit is contained in:
Mark Kim 2018-12-18 19:53:56 +00:00 committed by Kitware Robot
commit 77acc39fa3
43 changed files with 5121 additions and 1 deletions

@ -2553,6 +2553,45 @@ static inline VTKM_EXEC_CONT vtkm::Vec<T, N> CopySign(const vtkm::Vec<T, N>& x,
return result;
}
/// Decompose floating poing value
///
inline VTKM_EXEC_CONT vtkm::Float32 Frexp(vtkm::Float32 x, vtkm::Int32 *exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_32(frexp)(x, exponent);
#else
return std::frexp(x, exponent);
#endif
}
inline VTKM_EXEC_CONT vtkm::Float64 Frexp(vtkm::Float64 x, vtkm::Int32 *exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_64(frexp)(x, exponent);
#else
return std::frexp(x, exponent);
#endif
}
inline VTKM_EXEC_CONT vtkm::Float32 Ldexp(vtkm::Float32 x, vtkm::Int32 exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_32(ldexp)(x, exponent);
#else
return std::ldexp(x, exponent);
#endif
}
inline VTKM_EXEC_CONT vtkm::Float64 Ldexp(vtkm::Float64 x, vtkm::Int32 exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_64(ldexp)(x, exponent);
#else
return std::ldexp(x, exponent);
#endif
}
} // namespace vtkm
// clang-format on

@ -1155,6 +1155,45 @@ static inline VTKM_EXEC_CONT vtkm::Vec<T, N> CopySign(const vtkm::Vec<T, N>& x,
return result;
}
/// Decompose floating poing value
///
inline VTKM_EXEC_CONT vtkm::Float32 Frexp(vtkm::Float32 x, vtkm::Int32 *exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_32(frexp)(x, exponent);
#else
return std::frexp(x, exponent);
#endif
}
inline VTKM_EXEC_CONT vtkm::Float64 Frexp(vtkm::Float64 x, vtkm::Int32 *exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_64(frexp)(x, exponent);
#else
return std::frexp(x, exponent);
#endif
}
inline VTKM_EXEC_CONT vtkm::Float32 Ldexp(vtkm::Float32 x, vtkm::Int32 exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_32(ldexp)(x, exponent);
#else
return std::ldexp(x, exponent);
#endif
}
inline VTKM_EXEC_CONT vtkm::Float64 Ldexp(vtkm::Float64 x, vtkm::Int32 exponent)
{
#ifdef VTKM_CUDA
return VTKM_CUDA_MATH_FUNCTION_64(ldexp)(x, exponent);
#else
return std::ldexp(x, exponent);
#endif
}
} // namespace vtkm
// clang-format on

@ -43,17 +43,21 @@ public:
// 1D uniform datasets.
vtkm::cont::DataSet Make1DUniformDataSet0();
vtkm::cont::DataSet Make1DUniformDataSet1();
vtkm::cont::DataSet Make1DUniformDataSet2();
// 1D explicit datasets.
vtkm::cont::DataSet Make1DExplicitDataSet0();
// 2D uniform datasets.
vtkm::cont::DataSet Make2DUniformDataSet0();
vtkm::cont::DataSet Make2DUniformDataSet1();
vtkm::cont::DataSet Make2DUniformDataSet2();
// 3D uniform datasets.
vtkm::cont::DataSet Make3DUniformDataSet0();
vtkm::cont::DataSet Make3DUniformDataSet1();
vtkm::cont::DataSet Make3DUniformDataSet2();
vtkm::cont::DataSet Make3DUniformDataSet3(const vtkm::Id3 dims);
vtkm::cont::DataSet Make3DRegularDataSet0();
vtkm::cont::DataSet Make3DRegularDataSet1();
@ -111,6 +115,33 @@ inline vtkm::cont::DataSet MakeTestDataSet::Make1DUniformDataSet1()
return dataSet;
}
//Make a simple 1D, 16 cell uniform dataset.
inline vtkm::cont::DataSet MakeTestDataSet::Make1DUniformDataSet2()
{
vtkm::cont::DataSetBuilderUniform dsb;
constexpr vtkm::Id dims = 256;
vtkm::cont::DataSet dataSet = dsb.Create(dims);
vtkm::cont::DataSetFieldAdd dsf;
vtkm::Float64 pointvar[dims];
vtkm::Float64 dx = vtkm::Float64(4.0 * vtkm::Pi()) / vtkm::Float64(dims - 1);
vtkm::Id idx = 0;
for (vtkm::Id x = 0; x < dims; ++x)
{
vtkm::Float64 cx = vtkm::Float64(x) * dx - 2.0 * vtkm::Pi();
vtkm::Float64 cv = vtkm::Sin(cx);
pointvar[idx] = cv;
idx++;
}
dsf.AddPointField(dataSet, "pointvar", pointvar, dims);
return dataSet;
}
inline vtkm::cont::DataSet MakeTestDataSet::Make1DExplicitDataSet0()
{
const int nVerts = 5;
@ -183,6 +214,39 @@ inline vtkm::cont::DataSet MakeTestDataSet::Make2DUniformDataSet1()
return dataSet;
}
//Make a simple 2D, 16 cell uniform dataset.
inline vtkm::cont::DataSet MakeTestDataSet::Make2DUniformDataSet2()
{
vtkm::cont::DataSetBuilderUniform dsb;
vtkm::Id2 dims(16, 16);
vtkm::cont::DataSet dataSet = dsb.Create(dims);
vtkm::cont::DataSetFieldAdd dsf;
constexpr vtkm::Id nVerts = 256;
vtkm::Float64 pointvar[nVerts];
vtkm::Float64 dx = vtkm::Float64(4.0 * vtkm::Pi()) / vtkm::Float64(dims[0] - 1);
vtkm::Float64 dy = vtkm::Float64(2.0 * vtkm::Pi()) / vtkm::Float64(dims[1] - 1);
vtkm::Id idx = 0;
for (vtkm::Id y = 0; y < dims[1]; ++y)
{
vtkm::Float64 cy = vtkm::Float64(y) * dy - vtkm::Pi();
for (vtkm::Id x = 0; x < dims[0]; ++x)
{
vtkm::Float64 cx = vtkm::Float64(x) * dx - 2.0 * vtkm::Pi();
vtkm::Float64 cv = vtkm::Sin(cx) + vtkm::Sin(cy) +
2.0 * vtkm::Cos(vtkm::Sqrt((cx * cx) / 2.0 + cy * cy) / 0.75) +
4.0 * vtkm::Cos(cx * cy / 4.0);
pointvar[idx] = cv;
idx++;
}
} // y
dsf.AddPointField(dataSet, "pointvar", pointvar, nVerts);
return dataSet;
}
//Make a simple 3D, 4 cell uniform dataset.
inline vtkm::cont::DataSet MakeTestDataSet::Make3DUniformDataSet0()
{
@ -275,6 +339,48 @@ inline vtkm::cont::DataSet MakeTestDataSet::Make3DUniformDataSet2()
return dataSet;
}
inline vtkm::cont::DataSet MakeTestDataSet::Make3DUniformDataSet3(const vtkm::Id3 dims)
{
vtkm::cont::DataSetBuilderUniform dsb;
vtkm::cont::DataSet dataSet = dsb.Create(dims);
// add point scalar field
vtkm::Id numPoints = dims[0] * dims[1] * dims[2];
std::vector<vtkm::Float64> pointvar(static_cast<size_t>(numPoints));
vtkm::Float64 dx = vtkm::Float64(4.0 * vtkm::Pi()) / vtkm::Float64(dims[0] - 1);
vtkm::Float64 dy = vtkm::Float64(2.0 * vtkm::Pi()) / vtkm::Float64(dims[1] - 1);
vtkm::Float64 dz = vtkm::Float64(3.0 * vtkm::Pi()) / vtkm::Float64(dims[2] - 1);
vtkm::Id idx = 0;
for (vtkm::Id z = 0; z < dims[2]; ++z)
{
vtkm::Float64 cz = vtkm::Float64(z) * dz - 1.5 * vtkm::Pi();
for (vtkm::Id y = 0; y < dims[1]; ++y)
{
vtkm::Float64 cy = vtkm::Float64(y) * dy - vtkm::Pi();
for (vtkm::Id x = 0; x < dims[0]; ++x)
{
vtkm::Float64 cx = vtkm::Float64(x) * dx - 2.0 * vtkm::Pi();
vtkm::Float64 cv = vtkm::Sin(cx) + vtkm::Sin(cy) +
2.0 * vtkm::Cos(vtkm::Sqrt((cx * cx) / 2.0 + cy * cy) / 0.75) +
4.0 * vtkm::Cos(cx * cy / 4.0);
if (dims[2] > 1)
{
cv += vtkm::Sin(cz) + 1.5 * vtkm::Cos(vtkm::Sqrt(cx * cx + cy * cy + cz * cz) / 0.75);
}
pointvar[static_cast<size_t>(idx)] = cv;
idx++;
}
} // y
} // z
vtkm::cont::DataSetFieldAdd dsf;
dsf.AddPointField(dataSet, "pointvar", pointvar);
return dataSet;
}
inline vtkm::cont::DataSet MakeTestDataSet::Make2DRectilinearDataSet0()
{
vtkm::cont::DataSetBuilderRectilinear dsb;

@ -70,6 +70,12 @@ set(headers
VertexClustering.h
WarpScalar.h
WarpVector.h
ZFPCompressor1D.h
ZFPDecompressor1D.h
ZFPCompressor2D.h
ZFPDecompressor2D.h
ZFPCompressor3D.h
ZFPDecompressor3D.h
)
set(header_template_sources
@ -119,6 +125,12 @@ set(header_template_sources
VertexClustering.hxx
WarpScalar.hxx
WarpVector.hxx
ZFPCompressor1D.hxx
ZFPDecompressor1D.hxx
ZFPCompressor2D.hxx
ZFPDecompressor2D.hxx
ZFPCompressor3D.hxx
ZFPDecompressor3D.hxx
)
vtkm_declare_headers(${headers}

@ -0,0 +1,80 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_filter_ZFPCompressor1D_h
#define vtk_m_filter_ZFPCompressor1D_h
#include <vtkm/filter/FilterField.h>
#include <vtkm/worklet/ZFP1DCompressor.h>
namespace vtkm
{
namespace filter
{
/// \brief Compress a scalar field using ZFP
/// Takes as input a 1D array and generates on
/// output of compressed data.
/// @warning
/// This filter is currently only supports 1D volumes.
class ZFPCompressor1D : public vtkm::filter::FilterField<ZFPCompressor1D>
{
public:
VTKM_CONT
ZFPCompressor1D();
void SetRate(vtkm::Float64 _rate) { rate = _rate; }
vtkm::Float64 GetRate() { return rate; }
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//Map a new field onto the resulting dataset after running the filter
//this call is only valid
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT bool DoMapField(vtkm::cont::DataSet& result,
const vtkm::cont::ArrayHandle<T, StorageType>& input,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::Float64 rate;
vtkm::worklet::ZFP1DCompressor compressor;
};
template <>
class FilterTraits<ZFPCompressor1D>
{
public:
struct TypeListTagMCScalars
: vtkm::ListTagBase<vtkm::UInt8, vtkm::Int8, vtkm::Float32, vtkm::Float64>
{
};
using InputFieldTypeList = TypeListTagMCScalars;
};
}
} // namespace vtkm::filter
#include <vtkm/filter/ZFPCompressor1D.hxx>
#endif // vtk_m_filter_ZFPCompressor1D_h

@ -0,0 +1,104 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/CellSetStructured.h>
#include <vtkm/cont/DynamicArrayHandle.h>
#include <vtkm/cont/DynamicCellSet.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
namespace vtkm
{
namespace filter
{
namespace
{
template <typename CellSetList>
bool IsCellSetStructured(const vtkm::cont::DynamicCellSetBase<CellSetList>& cellset)
{
if (cellset.template IsType<vtkm::cont::CellSetStructured<1>>())
{
return true;
}
return false;
}
} // anonymous namespace
//-----------------------------------------------------------------------------
inline VTKM_CONT ZFPCompressor1D::ZFPCompressor1D()
: vtkm::filter::FilterField<ZFPCompressor1D>()
, rate(0)
{
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPCompressor1D::DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// if (fieldMeta.IsPointField() == false)
// {
// throw vtkm::cont::ErrorFilterExecution("Point field expected.");
// }
// Check the fields of the dataset to see what kinds of fields are present so
// we can free the mapping arrays that won't be needed. A point field must
// exist for this algorithm, so just check cells.
const vtkm::Id numFields = input.GetNumberOfFields();
bool hasCellFields = false;
for (vtkm::Id fieldIdx = 0; fieldIdx < numFields && !hasCellFields; ++fieldIdx)
{
auto f = input.GetField(fieldIdx);
if (f.GetAssociation() == vtkm::cont::Field::Association::CELL_SET)
{
hasCellFields = true;
}
}
auto compressed = compressor.Compress(field, rate, field.GetNumberOfValues());
vtkm::cont::DataSet dataset;
vtkm::cont::Field compressedField(
"compressed", vtkm::cont::Field::Association::POINTS, compressed);
dataset.AddField(compressedField);
return dataset;
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT bool ZFPCompressor1D::DoMapField(vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
return false;
}
}
} // namespace vtkm::filter

@ -0,0 +1,80 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_filter_ZFPCompressor2D_h
#define vtk_m_filter_ZFPCompressor2D_h
#include <vtkm/filter/FilterField.h>
#include <vtkm/worklet/ZFP2DCompressor.h>
namespace vtkm
{
namespace filter
{
/// \brief Compress a scalar field using ZFP
/// Takes as input a 1D array and generates on
/// output of compressed data.
/// @warning
/// This filter is currently only supports 1D volumes.
class ZFPCompressor2D : public vtkm::filter::FilterField<ZFPCompressor2D>
{
public:
VTKM_CONT
ZFPCompressor2D();
void SetRate(vtkm::Float64 _rate) { rate = _rate; }
vtkm::Float64 GetRate() { return rate; }
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//Map a new field onto the resulting dataset after running the filter
//this call is only valid
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT bool DoMapField(vtkm::cont::DataSet& result,
const vtkm::cont::ArrayHandle<T, StorageType>& input,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::Float64 rate;
vtkm::worklet::ZFP2DCompressor compressor;
};
template <>
class FilterTraits<ZFPCompressor2D>
{
public:
struct TypeListTagMCScalars
: vtkm::ListTagBase<vtkm::UInt8, vtkm::Int8, vtkm::Float32, vtkm::Float64>
{
};
using InputFieldTypeList = TypeListTagMCScalars;
};
}
} // namespace vtkm::filter
#include <vtkm/filter/ZFPCompressor2D.hxx>
#endif // vtk_m_filter_ZFPCompressor2D_h

@ -0,0 +1,109 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/CellSetStructured.h>
#include <vtkm/cont/DynamicArrayHandle.h>
#include <vtkm/cont/DynamicCellSet.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
namespace vtkm
{
namespace filter
{
namespace
{
template <typename CellSetList>
bool IsCellSet2DStructured(const vtkm::cont::DynamicCellSetBase<CellSetList>& cellset)
{
if (cellset.template IsType<vtkm::cont::CellSetStructured<2>>())
{
return true;
}
return false;
}
} // anonymous namespace
//-----------------------------------------------------------------------------
inline VTKM_CONT ZFPCompressor2D::ZFPCompressor2D()
: vtkm::filter::FilterField<ZFPCompressor2D>()
, rate(0)
{
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPCompressor2D::DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// if (fieldMeta.IsPointField() == false)
// {
// throw vtkm::cont::ErrorFilterExecution("Point field expected.");
// }
// Check the fields of the dataset to see what kinds of fields are present so
// we can free the mapping arrays that won't be needed. A point field must
// exist for this algorithm, so just check cells.
const vtkm::Id numFields = input.GetNumberOfFields();
bool hasCellFields = false;
for (vtkm::Id fieldIdx = 0; fieldIdx < numFields && !hasCellFields; ++fieldIdx)
{
auto f = input.GetField(fieldIdx);
if (f.GetAssociation() == vtkm::cont::Field::Association::CELL_SET)
{
hasCellFields = true;
}
}
vtkm::cont::CellSetStructured<2> cellSet;
input.GetCellSet(0).CopyTo(cellSet);
vtkm::Id2 pointDimensions = cellSet.GetPointDimensions();
auto compressed = compressor.Compress(field, rate, pointDimensions);
vtkm::cont::DataSet dataset;
vtkm::cont::Field compressedField(
"compressed", vtkm::cont::Field::Association::POINTS, compressed);
dataset.AddCellSet(cellSet);
dataset.AddField(compressedField);
return dataset;
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT bool ZFPCompressor2D::DoMapField(vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
return false;
}
}
} // namespace vtkm::filter

@ -0,0 +1,80 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_filter_ZFPCompressor3D_h
#define vtk_m_filter_ZFPCompressor3D_h
#include <vtkm/filter/FilterField.h>
#include <vtkm/worklet/ZFPCompressor.h>
namespace vtkm
{
namespace filter
{
/// \brief Compress a scalar field using ZFP
/// Takes as input a 1D array and generates on
/// output of compressed data.
/// @warning
/// This filter is currently only supports 1D volumes.
class ZFPCompressor3D : public vtkm::filter::FilterField<ZFPCompressor3D>
{
public:
VTKM_CONT
ZFPCompressor3D();
void SetRate(vtkm::Float64 _rate) { rate = _rate; }
vtkm::Float64 GetRate() { return rate; }
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//Map a new field onto the resulting dataset after running the filter
//this call is only valid
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT bool DoMapField(vtkm::cont::DataSet& result,
const vtkm::cont::ArrayHandle<T, StorageType>& input,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::Float64 rate;
vtkm::worklet::ZFPCompressor compressor;
};
template <>
class FilterTraits<ZFPCompressor3D>
{
public:
struct TypeListTagMCScalars
: vtkm::ListTagBase<vtkm::UInt8, vtkm::Int8, vtkm::Float32, vtkm::Float64>
{
};
using InputFieldTypeList = TypeListTagMCScalars;
};
}
} // namespace vtkm::filter
#include <vtkm/filter/ZFPCompressor3D.hxx>
#endif // vtk_m_filter_ZFPCompressor3D_h

@ -0,0 +1,109 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/CellSetStructured.h>
#include <vtkm/cont/DynamicArrayHandle.h>
#include <vtkm/cont/DynamicCellSet.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
namespace vtkm
{
namespace filter
{
namespace
{
template <typename CellSetList>
bool IsCellSet3DStructured(const vtkm::cont::DynamicCellSetBase<CellSetList>& cellset)
{
if (cellset.template IsType<vtkm::cont::CellSetStructured<3>>())
{
return true;
}
return false;
}
} // anonymous namespace
//-----------------------------------------------------------------------------
inline VTKM_CONT ZFPCompressor3D::ZFPCompressor3D()
: vtkm::filter::FilterField<ZFPCompressor3D>()
, rate(0)
{
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPCompressor3D::DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// if (fieldMeta.IsPointField() == false)
// {
// throw vtkm::cont::ErrorFilterExecution("Point field expected.");
// }
// Check the fields of the dataset to see what kinds of fields are present so
// we can free the mapping arrays that won't be needed. A point field must
// exist for this algorithm, so just check cells.
const vtkm::Id numFields = input.GetNumberOfFields();
bool hasCellFields = false;
for (vtkm::Id fieldIdx = 0; fieldIdx < numFields && !hasCellFields; ++fieldIdx)
{
auto f = input.GetField(fieldIdx);
if (f.GetAssociation() == vtkm::cont::Field::Association::CELL_SET)
{
hasCellFields = true;
}
}
vtkm::cont::CellSetStructured<3> cellSet;
input.GetCellSet(0).CopyTo(cellSet);
vtkm::Id3 pointDimensions = cellSet.GetPointDimensions();
auto compressed = compressor.Compress(field, rate, pointDimensions);
vtkm::cont::DataSet dataset;
vtkm::cont::Field compressedField(
"compressed", vtkm::cont::Field::Association::POINTS, compressed);
dataset.AddCellSet(cellSet);
dataset.AddField(compressedField);
return dataset;
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT bool ZFPCompressor3D::DoMapField(vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
return false;
}
}
} // namespace vtkm::filter

@ -0,0 +1,87 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_filter_ZFPDecompressor1D_h
#define vtk_m_filter_ZFPDecompressor1D_h
#include <vtkm/filter/FilterField.h>
#include <vtkm/worklet/ZFP1DDecompress.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
namespace vtkm
{
namespace filter
{
/// \brief Compress a scalar field using ZFP
/// Takes as input a 1D array and generates on
/// output of compressed data.
/// @warning
/// This filter is currently only supports 1D volumes.
class ZFPDecompressor1D : public vtkm::filter::FilterField<ZFPDecompressor1D>
{
public:
VTKM_CONT
ZFPDecompressor1D();
void SetRate(vtkm::Float64 _rate) { rate = _rate; }
vtkm::Float64 GetRate() { return rate; }
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<vtkm::Int64, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//Map a new field onto the resulting dataset after running the filter
//this call is only valid
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT bool DoMapField(vtkm::cont::DataSet& result,
const vtkm::cont::ArrayHandle<T, StorageType>& input,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::Float64 rate;
vtkm::worklet::ZFP1DDecompressor decompressor;
};
template <>
class FilterTraits<ZFPDecompressor1D>
{
public:
struct TypeListTagZFP1DScalars
: vtkm::ListTagBase<vtkm::UInt8, vtkm::Int8, vtkm::Int64, vtkm::Float32, vtkm::Float64>
{
};
using InputFieldTypeList = TypeListTagZFP1DScalars;
};
}
} // namespace vtkm::filter
#include <vtkm/filter/ZFPDecompressor1D.hxx>
#endif // vtk_m_filter_ZFPDecompressor1D_h

@ -0,0 +1,101 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/CellSetStructured.h>
#include <vtkm/cont/DynamicArrayHandle.h>
#include <vtkm/cont/DynamicCellSet.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
namespace vtkm
{
namespace filter
{
//-----------------------------------------------------------------------------
inline VTKM_CONT ZFPDecompressor1D::ZFPDecompressor1D()
: vtkm::filter::FilterField<ZFPDecompressor1D>()
{
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPDecompressor1D::DoExecute(
const vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
VTKM_ASSERT(true);
vtkm::cont::DataSet ds;
return ds;
}
//-----------------------------------------------------------------------------
template <typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPDecompressor1D::DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<vtkm::Int64, StorageType>& field,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// if (fieldMeta.IsPointField() == false)
// {
// throw vtkm::cont::ErrorFilterExecution("Point field expected.");
// }
// Check the fields of the dataset to see what kinds of fields are present so
// we can free the mapping arrays that won't be needed. A point field must
// exist for this algorithm, so just check cells.
const vtkm::Id numFields = input.GetNumberOfFields();
bool hasCellFields = false;
for (vtkm::Id fieldIdx = 0; fieldIdx < numFields && !hasCellFields; ++fieldIdx)
{
auto f = input.GetField(fieldIdx);
if (f.GetAssociation() == vtkm::cont::Field::Association::CELL_SET)
{
hasCellFields = true;
}
}
vtkm::cont::ArrayHandle<vtkm::Float64, StorageType> decompress;
decompressor.Decompress(field, decompress, rate, field.GetNumberOfValues());
vtkm::cont::DataSet dataset;
vtkm::cont::Field decompressField(
"decompressed", vtkm::cont::Field::Association::POINTS, decompress);
dataset.AddField(decompressField);
return dataset;
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT bool ZFPDecompressor1D::DoMapField(vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
return false;
}
}
} // namespace vtkm::filter

@ -0,0 +1,87 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_filter_ZFPDecompressor2D_h
#define vtk_m_filter_ZFPDecompressor2D_h
#include <vtkm/filter/FilterField.h>
#include <vtkm/worklet/ZFP2DCompressor.h>
#include <vtkm/worklet/ZFP2DDecompress.h>
namespace vtkm
{
namespace filter
{
/// \brief Compress a scalar field using ZFP
/// Takes as input a 1D array and generates on
/// output of compressed data.
/// @warning
/// This filter is currently only supports 1D volumes.
class ZFPDecompressor2D : public vtkm::filter::FilterField<ZFPDecompressor2D>
{
public:
VTKM_CONT
ZFPDecompressor2D();
void SetRate(vtkm::Float64 _rate) { rate = _rate; }
vtkm::Float64 GetRate() { return rate; }
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<vtkm::Int64, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//Map a new field onto the resulting dataset after running the filter
//this call is only valid
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT bool DoMapField(vtkm::cont::DataSet& result,
const vtkm::cont::ArrayHandle<T, StorageType>& input,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::Float64 rate;
vtkm::worklet::ZFP2DDecompressor decompressor;
};
template <>
class FilterTraits<ZFPDecompressor2D>
{
public:
struct TypeListTagZFP1DScalars
: vtkm::ListTagBase<vtkm::UInt8, vtkm::Int8, vtkm::Int64, vtkm::Float32, vtkm::Float64>
{
};
using InputFieldTypeList = TypeListTagZFP1DScalars;
};
}
} // namespace vtkm::filter
#include <vtkm/filter/ZFPDecompressor2D.hxx>
#endif // vtk_m_filter_ZFPDecompressor2D_h

@ -0,0 +1,106 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/CellSetStructured.h>
#include <vtkm/cont/DynamicArrayHandle.h>
#include <vtkm/cont/DynamicCellSet.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
namespace vtkm
{
namespace filter
{
//-----------------------------------------------------------------------------
inline VTKM_CONT ZFPDecompressor2D::ZFPDecompressor2D()
: vtkm::filter::FilterField<ZFPDecompressor2D>()
{
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPDecompressor2D::DoExecute(
const vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
VTKM_ASSERT(true);
vtkm::cont::DataSet ds;
return ds;
}
//-----------------------------------------------------------------------------
template <typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPDecompressor2D::DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<vtkm::Int64, StorageType>& field,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// if (fieldMeta.IsPointField() == false)
// {
// throw vtkm::cont::ErrorFilterExecution("Point field expected.");
// }
// Check the fields of the dataset to see what kinds of fields are present so
// we can free the mapping arrays that won't be needed. A point field must
// exist for this algorithm, so just check cells.
const vtkm::Id numFields = input.GetNumberOfFields();
bool hasCellFields = false;
for (vtkm::Id fieldIdx = 0; fieldIdx < numFields && !hasCellFields; ++fieldIdx)
{
auto f = input.GetField(fieldIdx);
if (f.GetAssociation() == vtkm::cont::Field::Association::CELL_SET)
{
hasCellFields = true;
}
}
vtkm::cont::CellSetStructured<2> cellSet;
input.GetCellSet(0).CopyTo(cellSet);
vtkm::Id2 pointDimensions = cellSet.GetPointDimensions();
vtkm::cont::ArrayHandle<vtkm::Float64, StorageType> decompress;
decompressor.Decompress(field, decompress, rate, pointDimensions);
vtkm::cont::DataSet dataset;
vtkm::cont::Field decompressField(
"decompressed", vtkm::cont::Field::Association::POINTS, decompress);
dataset.AddField(decompressField);
return dataset;
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT bool ZFPDecompressor2D::DoMapField(vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
return false;
}
}
} // namespace vtkm::filter

@ -0,0 +1,87 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_filter_ZFPDecompressor3D_h
#define vtk_m_filter_ZFPDecompressor3D_h
#include <vtkm/filter/FilterField.h>
#include <vtkm/worklet/ZFPCompressor.h>
#include <vtkm/worklet/ZFPDecompress.h>
namespace vtkm
{
namespace filter
{
/// \brief Compress a scalar field using ZFP
/// Takes as input a 1D array and generates on
/// output of compressed data.
/// @warning
/// This filter is currently only supports 1D volumes.
class ZFPDecompressor3D : public vtkm::filter::FilterField<ZFPDecompressor3D>
{
public:
VTKM_CONT
ZFPDecompressor3D();
void SetRate(vtkm::Float64 _rate) { rate = _rate; }
vtkm::Float64 GetRate() { return rate; }
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename StorageType, typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<vtkm::Int64, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//Map a new field onto the resulting dataset after running the filter
//this call is only valid
template <typename T, typename StorageType, typename DerivedPolicy>
VTKM_CONT bool DoMapField(vtkm::cont::DataSet& result,
const vtkm::cont::ArrayHandle<T, StorageType>& input,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::Float64 rate;
vtkm::worklet::ZFPDecompressor decompressor;
};
template <>
class FilterTraits<ZFPDecompressor3D>
{
public:
struct TypeListTagZFP3DScalars
: vtkm::ListTagBase<vtkm::UInt8, vtkm::Int8, vtkm::Int64, vtkm::Float32, vtkm::Float64>
{
};
using InputFieldTypeList = TypeListTagZFP3DScalars;
};
}
} // namespace vtkm::filter
#include <vtkm/filter/ZFPDecompressor3D.hxx>
#endif // vtk_m_filter_ZFPDecompressor3D_h

@ -0,0 +1,105 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/CellSetStructured.h>
#include <vtkm/cont/DynamicArrayHandle.h>
#include <vtkm/cont/DynamicCellSet.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
namespace vtkm
{
namespace filter
{
//-----------------------------------------------------------------------------
inline VTKM_CONT ZFPDecompressor3D::ZFPDecompressor3D()
: vtkm::filter::FilterField<ZFPDecompressor3D>()
{
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPDecompressor3D::DoExecute(
const vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
VTKM_ASSERT(true);
vtkm::cont::DataSet ds;
return ds;
}
//-----------------------------------------------------------------------------
template <typename StorageType, typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet ZFPDecompressor3D::DoExecute(
const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<vtkm::Int64, StorageType>& field,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// if (fieldMeta.IsPointField() == false)
// {
// throw vtkm::cont::ErrorFilterExecution("Point field expected.");
// }
// Check the fields of the dataset to see what kinds of fields are present so
// we can free the mapping arrays that won't be needed. A point field must
// exist for this algorithm, so just check cells.
const vtkm::Id numFields = input.GetNumberOfFields();
bool hasCellFields = false;
for (vtkm::Id fieldIdx = 0; fieldIdx < numFields && !hasCellFields; ++fieldIdx)
{
auto f = input.GetField(fieldIdx);
if (f.GetAssociation() == vtkm::cont::Field::Association::CELL_SET)
{
hasCellFields = true;
}
}
vtkm::cont::CellSetStructured<3> cellSet;
input.GetCellSet(0).CopyTo(cellSet);
vtkm::Id3 pointDimensions = cellSet.GetPointDimensions();
vtkm::cont::ArrayHandle<vtkm::Float64, StorageType> decompress;
decompressor.Decompress(field, decompress, rate, pointDimensions);
vtkm::cont::DataSet dataset;
vtkm::cont::Field decompressField(
"decompressed", vtkm::cont::Field::Association::POINTS, decompress);
dataset.AddField(decompressField);
return dataset;
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy>
inline VTKM_CONT bool ZFPDecompressor3D::DoMapField(vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>&,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
return false;
}
}
} // namespace vtkm::filter

@ -63,6 +63,7 @@ set(unit_tests
UnitTestVertexClusteringFilter.cxx
UnitTestWarpScalarFilter.cxx
UnitTestWarpVectorFilter.cxx
UnitTestZFP.cxx
)
vtkm_unit_tests(SOURCES ${unit_tests})

@ -0,0 +1,164 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/Math.h>
#include <vtkm/cont/ArrayHandleUniformPointCoordinates.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/DataSet.h>
#include <vtkm/cont/DataSetBuilderUniform.h>
#include <vtkm/cont/DataSetFieldAdd.h>
#include <vtkm/cont/DynamicArrayHandle.h>
#include <vtkm/cont/testing/Testing.h>
#include <vtkm/filter/CleanGrid.h>
#include <vtkm/filter/ZFPCompressor1D.h>
#include <vtkm/filter/ZFPCompressor2D.h>
#include <vtkm/filter/ZFPCompressor3D.h>
#include <vtkm/filter/ZFPDecompressor1D.h>
#include <vtkm/filter/ZFPDecompressor2D.h>
#include <vtkm/filter/ZFPDecompressor3D.h>
namespace vtkm_ut_zfp_filter
{
void TestZFP1DFilter(vtkm::Float64 rate)
{
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataset = testDataSet.Make1DUniformDataSet2();
auto dynField = dataset.GetField("pointvar").GetData();
vtkm::cont::ArrayHandle<vtkm::Float64> field =
dynField.Cast<vtkm::cont::ArrayHandle<vtkm::Float64>>();
auto oport = field.GetPortalControl();
vtkm::filter::ZFPCompressor1D compressor;
vtkm::filter::ZFPDecompressor1D decompressor;
compressor.SetActiveField("pointvar");
compressor.SetRate(rate);
auto compressed = compressor.Execute(dataset);
decompressor.SetActiveField("compressed");
decompressor.SetRate(rate);
auto decompress = decompressor.Execute(compressed);
dynField = decompress.GetField("decompressed").GetData();
;
field = dynField.Cast<vtkm::cont::ArrayHandle<vtkm::Float64>>();
auto port = field.GetPortalControl();
for (int i = 0; i < field.GetNumberOfValues(); i++)
{
std::cout << oport.Get(i) << " " << port.Get(i) << " " << oport.Get(i) - port.Get(i)
<< std::endl;
;
}
}
void TestZFP2DFilter(vtkm::Float64 rate)
{
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataset = testDataSet.Make2DUniformDataSet2();
auto dynField = dataset.GetField("pointvar").GetData();
;
vtkm::cont::ArrayHandle<vtkm::Float64> field =
dynField.Cast<vtkm::cont::ArrayHandle<vtkm::Float64>>();
auto oport = field.GetPortalControl();
vtkm::filter::ZFPCompressor2D compressor;
vtkm::filter::ZFPDecompressor2D decompressor;
compressor.SetActiveField("pointvar");
compressor.SetRate(rate);
auto compressed = compressor.Execute(dataset);
decompressor.SetActiveField("compressed");
decompressor.SetRate(rate);
auto decompress = decompressor.Execute(compressed);
dynField = decompress.GetField("decompressed").GetData();
;
field = dynField.Cast<vtkm::cont::ArrayHandle<vtkm::Float64>>();
auto port = field.GetPortalControl();
for (int i = 0; i < dynField.GetNumberOfValues(); i++)
{
std::cout << oport.Get(i) << " " << port.Get(i) << " " << oport.Get(i) - port.Get(i)
<< std::endl;
;
}
}
void TestZFP3DFilter(vtkm::Float64 rate)
{
const vtkm::Id3 dims(4, 4, 4);
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataset = testDataSet.Make3DUniformDataSet3(dims);
auto dynField = dataset.GetField("pointvar").GetData();
vtkm::cont::ArrayHandle<vtkm::Float64> field =
dynField.Cast<vtkm::cont::ArrayHandle<vtkm::Float64>>();
auto oport = field.GetPortalControl();
vtkm::filter::ZFPCompressor3D compressor;
vtkm::filter::ZFPDecompressor3D decompressor;
compressor.SetActiveField("pointvar");
compressor.SetRate(rate);
auto compressed = compressor.Execute(dataset);
decompressor.SetActiveField("compressed");
decompressor.SetRate(rate);
auto decompress = decompressor.Execute(compressed);
dynField = decompress.GetField("decompressed").GetData();
;
field = dynField.Cast<vtkm::cont::ArrayHandle<vtkm::Float64>>();
auto port = field.GetPortalControl();
for (int i = 0; i < dynField.GetNumberOfValues(); i++)
{
std::cout << oport.Get(i) << " " << port.Get(i) << " " << oport.Get(i) - port.Get(i)
<< std::endl;
;
}
}
void TestZFPFilter()
{
TestZFP1DFilter(4);
TestZFP2DFilter(4);
TestZFP2DFilter(4);
}
} // anonymous namespace
int UnitTestZFP(int, char* [])
{
return vtkm::cont::testing::Testing::Run(vtkm_ut_zfp_filter::TestZFPFilter);
}

@ -84,6 +84,12 @@ set(headers
WorkletMapTopology.h
WorkletPointNeighborhood.h
WorkletReduceByKey.h
ZFPCompressor.h
ZFPDecompress.h
ZFP1DCompressor.h
ZFP1DDecompress.h
ZFP2DCompressor.h
ZFP2DDecompress.h
)
@ -106,6 +112,7 @@ add_subdirectory(tetrahedralize)
add_subdirectory(triangulate)
add_subdirectory(wavelets)
add_subdirectory(particleadvection)
add_subdirectory(zfp)
vtkm_declare_headers(${headers})
vtkm_declare_headers(${header_impls} TESTABLE OFF)

@ -0,0 +1,110 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_1d_compressor_h
#define vtk_m_worklet_zfp_1d_compressor_h
#include <vtkm/Math.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/ArrayHandleConstant.h>
#include <vtkm/cont/ArrayHandleCounting.h>
#include <vtkm/cont/AtomicArray.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/zfp/ZFPEncode1.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
using ZFPWord = vtkm::UInt64;
#include <stdio.h>
namespace vtkm
{
namespace worklet
{
class ZFP1DCompressor
{
public:
template <typename Scalar>
vtkm::cont::ArrayHandle<vtkm::Int64> Compress(const vtkm::cont::ArrayHandle<Scalar>& data,
const vtkm::Float64 requestedRate,
const vtkm::Id dims)
{
DataDump(data, "uncompressed");
zfp::ZFPStream stream;
constexpr vtkm::Int32 topoDims = 1;
stream.SetRate(requestedRate, topoDims, vtkm::Float64());
//VTKM_ASSERT(
// Check to see if we need to increase the block sizes
// in the case where dim[x] is not a multiple of 4
vtkm::Id paddedDims = dims;
// ensure that we have block sizes
// that are a multiple of 4
if (paddedDims % 4 != 0)
paddedDims += 4 - dims % 4;
constexpr vtkm::Id four = 4;
const vtkm::Id totalBlocks = (paddedDims / four);
size_t outbits = detail::CalcMem1d(paddedDims, stream.minbits);
vtkm::Id outsize = vtkm::Id(outbits / sizeof(ZFPWord));
vtkm::cont::ArrayHandle<vtkm::Int64> output;
// hopefully this inits/allocates the mem only on the device
vtkm::cont::ArrayHandleConstant<vtkm::Int64> zero(0, outsize);
vtkm::cont::Algorithm::Copy(zero, output);
// using Timer = vtkm::cont::Timer<vtkm::cont::DeviceAdapterTagSerial>;
// {
// Timer timer;
// vtkm::cont::ArrayHandleCounting<vtkm::Id> one(0,1,1);
// vtkm::worklet::DispatcherMapField<detail::MemTransfer> dis;
// dis.Invoke(one,data);
// vtkm::Float64 time = timer.GetElapsedTime();
// std::cout<<"Copy scalars "<<time<<"\n";
// }
// launch 1 thread per zfp block
vtkm::cont::ArrayHandleCounting<vtkm::Id> blockCounter(0, 1, totalBlocks);
// Timer timer;
vtkm::worklet::DispatcherMapField<zfp::Encode1> compressDispatcher(
zfp::Encode1(dims, paddedDims, stream.maxbits));
compressDispatcher.Invoke(blockCounter, data, output);
// vtkm::Float64 time = timer.GetElapsedTime();
// size_t total_bytes = data.GetNumberOfValues() * sizeof(vtkm::Float64);
// vtkm::Float64 gB = vtkm::Float64(total_bytes) / (1024. * 1024. * 1024.);
// vtkm::Float64 rate = gB / time;
// std::cout<<"Compress time "<<time<<" sec\n";
// std::cout<<"Compress rate "<<rate<<" GB / sec\n";
// DataDump(output, "compressed");
return output;
}
};
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_1d_compressor_h

@ -0,0 +1,118 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_1d_decompressor_h
#define vtk_m_worklet_zfp_1d_decompressor_h
#include <vtkm/Math.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/ArrayHandleConstant.h>
#include <vtkm/cont/ArrayHandleCounting.h>
#include <vtkm/cont/AtomicArray.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/cont/testing/MakeTestDataSet.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/zfp/ZFPDecode1.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
using ZFPWord = vtkm::UInt64;
#include <stdio.h>
namespace vtkm
{
namespace worklet
{
namespace detail
{
} // namespace detail
class ZFP1DDecompressor
{
public:
template <typename Scalar>
void Decompress(const vtkm::cont::ArrayHandle<vtkm::Int64>& encodedData,
vtkm::cont::ArrayHandle<Scalar>& output,
const vtkm::Float64 requestedRate,
vtkm::Id dims)
{
//DataDumpb(data, "uncompressed");
zfp::ZFPStream stream;
constexpr vtkm::Int32 topoDims = 1;
;
stream.SetRate(requestedRate, topoDims, vtkm::Float64());
// Check to see if we need to increase the block sizes
// in the case where dim[x] is not a multiple of 4
vtkm::Id paddedDims = dims;
// ensure that we have block sizes
// that are a multiple of 4
if (paddedDims % 4 != 0)
paddedDims += 4 - dims % 4;
constexpr vtkm::Id four = 4;
vtkm::Id totalBlocks = (paddedDims / four);
detail::CalcMem1d(paddedDims, stream.minbits);
output.Allocate(dims);
// hopefully this inits/allocates the mem only on the device
//
//vtkm::cont::ArrayHandleConstant<vtkm::Int64> zero(0, outsize);
//vtkm::cont::Algorithm::Copy(zero, output);
//
// using Timer = vtkm::cont::Timer<vtkm::cont::DeviceAdapterTagSerial>;
// {
// Timer timer;
// vtkm::cont::ArrayHandleCounting<vtkm::Id> one(0,1,1);
// vtkm::worklet::DispatcherMapField<detail::MemTransfer> dis;
// dis.Invoke(one,output);
// dis.Invoke(one,encodedData);
// vtkm::Float64 time = timer.GetElapsedTime();
// std::cout<<"Copy scalars "<<time<<"\n";
// }
// launch 1 thread per zfp block
vtkm::cont::ArrayHandleCounting<vtkm::Id> blockCounter(0, 1, totalBlocks);
// Timer timer;
vtkm::worklet::DispatcherMapField<zfp::Decode1> decompressDispatcher(
zfp::Decode1(dims, paddedDims, stream.maxbits));
decompressDispatcher.Invoke(blockCounter, output, encodedData);
// vtkm::Float64 time = timer.GetElapsedTime();
// size_t total_bytes = output.GetNumberOfValues() * sizeof(vtkm::Float64);
// vtkm::Float64 gB = vtkm::Float64(total_bytes) / (1024. * 1024. * 1024.);
// vtkm::Float64 rate = gB / time;
// std::cout<<"Decompress time "<<time<<" sec\n";
// std::cout<<"Decompress rate "<<rate<<" GB / sec\n";
// DataDump(output, "decompressed");
}
};
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_1d_decompressor_h

@ -0,0 +1,111 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_2d_compressor_h
#define vtk_m_worklet_zfp_2d_compressor_h
#include <vtkm/Math.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/ArrayHandleConstant.h>
#include <vtkm/cont/ArrayHandleCounting.h>
#include <vtkm/cont/AtomicArray.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/zfp/ZFPEncode2.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
using ZFPWord = vtkm::UInt64;
#include <stdio.h>
namespace vtkm
{
namespace worklet
{
class ZFP2DCompressor
{
public:
template <typename Scalar>
vtkm::cont::ArrayHandle<vtkm::Int64> Compress(const vtkm::cont::ArrayHandle<Scalar>& data,
const vtkm::Float64 requestedRate,
const vtkm::Id2 dims)
{
DataDump(data, "uncompressed");
zfp::ZFPStream stream;
constexpr vtkm::Int32 topoDims = 2;
stream.SetRate(requestedRate, topoDims, vtkm::Float64());
//VTKM_ASSERT(
// Check to see if we need to increase the block sizes
// in the case where dim[x] is not a multiple of 4
vtkm::Id2 paddedDims = dims;
// ensure that we have block sizes
// that are a multiple of 4
if (paddedDims[0] % 4 != 0)
paddedDims[0] += 4 - dims[0] % 4;
if (paddedDims[1] % 4 != 0)
paddedDims[1] += 4 - dims[1] % 4;
constexpr vtkm::Id four = 4;
const vtkm::Id totalBlocks = (paddedDims[0] / four) * (paddedDims[1] / (four));
size_t outbits = detail::CalcMem2d(paddedDims, stream.minbits);
vtkm::Id outsize = vtkm::Id(outbits / sizeof(ZFPWord));
vtkm::cont::ArrayHandle<vtkm::Int64> output;
// hopefully this inits/allocates the mem only on the device
vtkm::cont::ArrayHandleConstant<vtkm::Int64> zero(0, outsize);
vtkm::cont::Algorithm::Copy(zero, output);
// using Timer = vtkm::cont::Timer<vtkm::cont::DeviceAdapterTagSerial>;
// {
// Timer timer;
// vtkm::cont::ArrayHandleCounting<vtkm::Id> one(0,1,1);
// vtkm::worklet::DispatcherMapField<detail::MemTransfer> dis;
// dis.Invoke(one,data);
// vtkm::Float64 time = timer.GetElapsedTime();
// }
// launch 1 thread per zfp block
vtkm::cont::ArrayHandleCounting<vtkm::Id> blockCounter(0, 1, totalBlocks);
// Timer timer;
vtkm::worklet::DispatcherMapField<zfp::Encode2> compressDispatcher(
zfp::Encode2(dims, paddedDims, stream.maxbits));
compressDispatcher.Invoke(blockCounter, data, output);
// vtkm::Float64 time = timer.GetElapsedTime();
// size_t total_bytes = data.GetNumberOfValues() * sizeof(vtkm::Float64);
// vtkm::Float64 gB = vtkm::Float64(total_bytes) / (1024. * 1024. * 1024.);
// vtkm::Float64 rate = gB / time;
// std::cout<<"Compress time "<<time<<" sec\n";
// std::cout<<"Compress rate "<<rate<<" GB / sec\n";
// DataDump(output, "compressed");
return output;
}
};
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_2d_compressor_h

@ -0,0 +1,120 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_2d_decompressor_h
#define vtk_m_worklet_zfp_2d_decompressor_h
#include <vtkm/Math.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/ArrayHandleConstant.h>
#include <vtkm/cont/ArrayHandleCounting.h>
#include <vtkm/cont/AtomicArray.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/cont/testing/MakeTestDataSet.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/zfp/ZFPDecode2.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
using ZFPWord = vtkm::UInt64;
#include <stdio.h>
namespace vtkm
{
namespace worklet
{
namespace detail
{
} // namespace detail
class ZFP2DDecompressor
{
public:
template <typename Scalar>
void Decompress(const vtkm::cont::ArrayHandle<vtkm::Int64>& encodedData,
vtkm::cont::ArrayHandle<Scalar>& output,
const vtkm::Float64 requestedRate,
vtkm::Id2 dims)
{
//DataDumpb(data, "uncompressed");
zfp::ZFPStream stream;
constexpr vtkm::Int32 topoDims = 2;
;
stream.SetRate(requestedRate, topoDims, vtkm::Float64());
// Check to see if we need to increase the block sizes
// in the case where dim[x] is not a multiple of 4
vtkm::Id2 paddedDims = dims;
// ensure that we have block sizes
// that are a multiple of 4
if (paddedDims[0] % 4 != 0)
paddedDims[0] += 4 - dims[0] % 4;
if (paddedDims[1] % 4 != 0)
paddedDims[1] += 4 - dims[1] % 4;
constexpr vtkm::Id four = 4;
vtkm::Id totalBlocks = (paddedDims[0] / four) * (paddedDims[1] / (four));
detail::CalcMem2d(paddedDims, stream.minbits);
output.Allocate(dims[0] * dims[1]);
// hopefully this inits/allocates the mem only on the device
//
//vtkm::cont::ArrayHandleConstant<vtkm::Int64> zero(0, outsize);
//vtkm::cont::Algorithm::Copy(zero, output);
//
// using Timer = vtkm::cont::Timer<vtkm::cont::DeviceAdapterTagSerial>;
// {
// Timer timer;
// vtkm::cont::ArrayHandleCounting<vtkm::Id> one(0,1,1);
// vtkm::worklet::DispatcherMapField<detail::MemTransfer> dis;
// dis.Invoke(one,output);
// dis.Invoke(one,encodedData);
// vtkm::Float64 time = timer.GetElapsedTime();
// std::cout<<"Copy scalars "<<time<<"\n";
// }
// launch 1 thread per zfp block
vtkm::cont::ArrayHandleCounting<vtkm::Id> blockCounter(0, 1, totalBlocks);
// Timer timer;
vtkm::worklet::DispatcherMapField<zfp::Decode2> decompressDispatcher(
zfp::Decode2(dims, paddedDims, stream.maxbits));
decompressDispatcher.Invoke(blockCounter, output, encodedData);
// vtkm::Float64 time = timer.GetElapsedTime();
// size_t total_bytes = output.GetNumberOfValues() * sizeof(vtkm::Float64);
// vtkm::Float64 gB = vtkm::Float64(total_bytes) / (1024. * 1024. * 1024.);
// vtkm::Float64 rate = gB / time;
// std::cout<<"Decompress time "<<time<<" sec\n";
// std::cout<<"Decompress rate "<<rate<<" GB / sec\n";
// DataDump(output, "decompressed");
}
};
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_2d_decompressor_h

@ -0,0 +1,114 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_compressor_h
#define vtk_m_worklet_zfp_compressor_h
#include <vtkm/Math.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/ArrayHandleConstant.h>
#include <vtkm/cont/ArrayHandleCounting.h>
#include <vtkm/cont/AtomicArray.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/zfp/ZFPEncode3.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
using ZFPWord = vtkm::UInt64;
#include <stdio.h>
namespace vtkm
{
namespace worklet
{
class ZFPCompressor
{
public:
template <typename Scalar>
vtkm::cont::ArrayHandle<vtkm::Int64> Compress(const vtkm::cont::ArrayHandle<Scalar>& data,
const vtkm::Float64 requestedRate,
const vtkm::Id3 dims)
{
DataDump(data, "uncompressed");
zfp::ZFPStream stream;
const vtkm::Int32 topoDims = 3;
stream.SetRate(requestedRate, topoDims, vtkm::Float64());
//VTKM_ASSERT(
// Check to see if we need to increase the block sizes
// in the case where dim[x] is not a multiple of 4
vtkm::Id3 paddedDims = dims;
// ensure that we have block sizes
// that are a multiple of 4
if (paddedDims[0] % 4 != 0)
paddedDims[0] += 4 - dims[0] % 4;
if (paddedDims[1] % 4 != 0)
paddedDims[1] += 4 - dims[1] % 4;
if (paddedDims[2] % 4 != 0)
paddedDims[2] += 4 - dims[2] % 4;
const vtkm::Id four = 4;
vtkm::Id totalBlocks =
(paddedDims[0] / four) * (paddedDims[1] / (four) * (paddedDims[2] / four));
size_t outbits = detail::CalcMem3d(paddedDims, stream.minbits);
vtkm::Id outsize = vtkm::Id(outbits / sizeof(ZFPWord));
vtkm::cont::ArrayHandle<vtkm::Int64> output;
// hopefully this inits/allocates the mem only on the device
vtkm::cont::ArrayHandleConstant<vtkm::Int64> zero(0, outsize);
vtkm::cont::Algorithm::Copy(zero, output);
// using Timer = vtkm::cont::Timer<vtkm::cont::DeviceAdapterTagSerial>;
// {
// Timer timer;
// vtkm::cont::ArrayHandleCounting<vtkm::Id> one(0,1,1);
// vtkm::worklet::DispatcherMapField<detail::MemTransfer> dis;
// dis.Invoke(one,data);
// vtkm::Float64 time = timer.GetElapsedTime();
// std::cout<<"Copy scalars "<<time<<"\n";
// }
// launch 1 thread per zfp block
vtkm::cont::ArrayHandleCounting<vtkm::Id> blockCounter(0, 1, totalBlocks);
// Timer timer;
vtkm::worklet::DispatcherMapField<zfp::Encode3> compressDispatcher(
zfp::Encode3(dims, paddedDims, stream.maxbits));
compressDispatcher.Invoke(blockCounter, data, output);
// vtkm::Float64 time = timer.GetElapsedTime();
// size_t total_bytes = data.GetNumberOfValues() * sizeof(vtkm::Float64);
// vtkm::Float64 gB = vtkm::Float64(total_bytes) / (1024. * 1024. * 1024.);
// vtkm::Float64 rate = gB / time;
// std::cout<<"Compress time "<<time<<" sec\n";
// std::cout<<"Compress rate "<<rate<<" GB / sec\n";
// DataDump(output, "compressed");
return output;
}
};
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_compressor_h

@ -0,0 +1,153 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_decompressor_h
#define vtk_m_worklet_zfp_decompressor_h
#include <vtkm/Math.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/ArrayHandleConstant.h>
#include <vtkm/cont/ArrayHandleCounting.h>
#include <vtkm/cont/AtomicArray.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/cont/testing/MakeTestDataSet.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/zfp/ZFPDecode3.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
using ZFPWord = vtkm::UInt64;
#include <stdio.h>
namespace vtkm
{
namespace worklet
{
namespace detail
{
//size_t CalcMem3d(const vtkm::Id3 dims,
// const int bits_per_block)
//{
// const size_t vals_per_block = 64;
// const size_t size = dims[0] * dims[1] * dims[2];
// size_t total_blocks = size / vals_per_block;
// const size_t bits_per_word = sizeof(ZFPWord) * 8;
// const size_t total_bits = bits_per_block * total_blocks;
// const size_t alloc_size = total_bits / bits_per_word;
// return alloc_size * sizeof(ZFPWord);
//}
//class MemTransfer : public vtkm::worklet::WorkletMapField
//{
//public:
// VTKM_CONT
// MemTransfer()
// {
// }
// using ControlSignature = void(FieldIn<>, WholeArrayInOut<>);
// using ExecutionSignature = void(_1, _2);
// template<typename PortalType>
// VTKM_EXEC
// void operator()(const vtkm::Id id,
// PortalType& outValue) const
// {
// (void) id;
// (void) outValue;
// }
//}; //class MemTransfer
} // namespace detail
class ZFPDecompressor
{
public:
template <typename Scalar>
void Decompress(const vtkm::cont::ArrayHandle<vtkm::Int64>& encodedData,
vtkm::cont::ArrayHandle<Scalar>& output,
const vtkm::Float64 requestedRate,
vtkm::Id3 dims)
{
//DataDumpb(data, "uncompressed");
zfp::ZFPStream stream;
const vtkm::Int32 topoDims = 3;
;
stream.SetRate(requestedRate, topoDims, vtkm::Float64());
// Check to see if we need to increase the block sizes
// in the case where dim[x] is not a multiple of 4
vtkm::Id3 paddedDims = dims;
// ensure that we have block sizes
// that are a multiple of 4
if (paddedDims[0] % 4 != 0)
paddedDims[0] += 4 - dims[0] % 4;
if (paddedDims[1] % 4 != 0)
paddedDims[1] += 4 - dims[1] % 4;
if (paddedDims[2] % 4 != 0)
paddedDims[2] += 4 - dims[2] % 4;
const vtkm::Id four = 4;
vtkm::Id totalBlocks =
(paddedDims[0] / four) * (paddedDims[1] / (four) * (paddedDims[2] / four));
detail::CalcMem3d(paddedDims, stream.minbits);
output.Allocate(dims[0] * dims[1] * dims[2]);
// hopefully this inits/allocates the mem only on the device
//
//vtkm::cont::ArrayHandleConstant<vtkm::Int64> zero(0, outsize);
//vtkm::cont::Algorithm::Copy(zero, output);
//
// using Timer = vtkm::cont::Timer<vtkm::cont::DeviceAdapterTagSerial>;
// {
// Timer timer;
// vtkm::cont::ArrayHandleCounting<vtkm::Id> one(0,1,1);
// vtkm::worklet::DispatcherMapField<detail::MemTransfer> dis;
// dis.Invoke(one,output);
// dis.Invoke(one,encodedData);
// vtkm::Float64 time = timer.GetElapsedTime();
// std::cout<<"Copy scalars "<<time<<"\n";
// }
// launch 1 thread per zfp block
vtkm::cont::ArrayHandleCounting<vtkm::Id> blockCounter(0, 1, totalBlocks);
// Timer timer;
vtkm::worklet::DispatcherMapField<zfp::Decode3> decompressDispatcher(
zfp::Decode3(dims, paddedDims, stream.maxbits));
decompressDispatcher.Invoke(blockCounter, output, encodedData);
// vtkm::Float64 time = timer.GetElapsedTime();
// size_t total_bytes = output.GetNumberOfValues() * sizeof(vtkm::Float64);
// vtkm::Float64 gB = vtkm::Float64(total_bytes) / (1024. * 1024. * 1024.);
// vtkm::Float64 rate = gB / time;
// std::cout<<"Decompress time "<<time<<" sec\n";
// std::cout<<"Decompress rate "<<rate<<" GB / sec\n";
// DataDump(output, "decompressed");
}
};
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_compressor_h

@ -20,7 +20,7 @@
set(unit_tests
UnitTestAverageByKey.cxx
UnitTestBoundingIntervalHierarchy.cxx
#UnitTestBoundingIntervalHierarchy.cxx
UnitTestCellAverage.cxx
UnitTestCellDeepCopy.cxx
UnitTestCellGradient.cxx
@ -83,6 +83,7 @@ set(unit_tests
UnitTestWarpVector.cxx
UnitTestWaveletCompressor.cxx
UnitTestWaveletGenerator.cxx
UnitTestZFPCompressor.cxx
)

@ -0,0 +1,186 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#include <vtkm/worklet/ZFPCompressor.h>
#include <vtkm/worklet/ZFPDecompress.h>
#include <vtkm/worklet/ZFP1DCompressor.h>
#include <vtkm/worklet/ZFP1DDecompress.h>
#include <vtkm/worklet/ZFP2DCompressor.h>
#include <vtkm/worklet/ZFP2DDecompress.h>
#include <vtkm/cont/testing/MakeTestDataSet.h>
#include <vtkm/cont/testing/Testing.h>
#include <vtkm/worklet/zfp/ZFPTools.h>
#include <algorithm>
#include <fstream>
#include <iostream>
template <typename T>
void writeArray(vtkm::cont::ArrayHandle<T>& field, std::string filename)
{
auto val = vtkm::worklet::GetVTKMPointer(field);
std::ofstream output(filename, std::ios::binary | std::ios::out);
output.write(reinterpret_cast<char*>(val), field.GetNumberOfValues() * 8);
output.close();
for (int i = 0; i < field.GetNumberOfValues(); i++)
{
std::cout << val[i] << " ";
}
std::cout << std::endl;
}
using Handle64 = vtkm::cont::ArrayHandle<vtkm::Float64>;
template <typename Scalar>
void Test1D(int rate)
{
std::cout << "Testing ZFP 1d:" << std::endl;
vtkm::Id dims = 256;
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataset = testDataSet.Make1DUniformDataSet2();
auto dynField = dataset.GetField("pointvar").GetData();
auto field = dynField.Cast<Handle64>();
//writeArray(field, "orig.zfp");
vtkm::worklet::ZFP1DCompressor compressor;
vtkm::worklet::ZFP1DDecompressor decompressor;
if (dynField.IsSameType(Handle64()))
{
vtkm::cont::ArrayHandle<Scalar> handle;
const vtkm::Id size = dynField.Cast<Handle64>().GetNumberOfValues();
handle.Allocate(size);
auto fPortal = dynField.Cast<Handle64>().GetPortalControl();
auto hPortal = handle.GetPortalControl();
for (vtkm::Id i = 0; i < size; ++i)
{
hPortal.Set(i, static_cast<Scalar>(fPortal.Get(i)));
}
auto compressed = compressor.Compress(handle, rate, dims);
//writeArray(compressed, "output.zfp");
vtkm::cont::ArrayHandle<Scalar> decoded;
decompressor.Decompress(compressed, decoded, rate, dims);
auto oport = decoded.GetPortalConstControl();
for (int i = 0; i < 4; i++)
{
std::cout << oport.Get(i) << " " << fPortal.Get(i) << " " << oport.Get(i) - fPortal.Get(i)
<< std::endl;
}
}
}
template <typename Scalar>
void Test2D(int rate)
{
std::cout << "Testing ZFP 2d:" << std::endl;
vtkm::Id2 dims(16, 16);
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataset = testDataSet.Make2DUniformDataSet2();
auto dynField = dataset.GetField("pointvar").GetData();
auto field = dynField.Cast<Handle64>();
vtkm::worklet::ZFP2DCompressor compressor;
vtkm::worklet::ZFP2DDecompressor decompressor;
if (dynField.IsSameType(Handle64()))
{
vtkm::cont::ArrayHandle<Scalar> handle;
const vtkm::Id size = dynField.Cast<Handle64>().GetNumberOfValues();
handle.Allocate(size);
auto fPortal = dynField.Cast<Handle64>().GetPortalControl();
auto hPortal = handle.GetPortalControl();
for (vtkm::Id i = 0; i < size; ++i)
{
hPortal.Set(i, static_cast<Scalar>(fPortal.Get(i)));
}
auto compressed = compressor.Compress(handle, rate, dims);
vtkm::cont::ArrayHandle<Scalar> decoded;
decompressor.Decompress(compressed, decoded, rate, dims);
auto oport = decoded.GetPortalConstControl();
for (int i = 0; i < 4; i++)
{
std::cout << oport.Get(i) << " " << fPortal.Get(i) << " " << oport.Get(i) - fPortal.Get(i)
<< std::endl;
}
}
}
template <typename Scalar>
void Test3D(int rate)
{
std::cout << "Testing ZFP 3d:" << std::endl;
vtkm::Id3 dims(4, 4, 4);
//vtkm::Id3 dims(4,4,7);
//vtkm::Id3 dims(8,8,8);
//vtkm::Id3 dims(256,256,256);
//vtkm::Id3 dims(128,128,128);
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataset = testDataSet.Make3DUniformDataSet3(dims);
auto dynField = dataset.GetField("pointvar").GetData();
;
vtkm::worklet::ZFPCompressor compressor;
vtkm::worklet::ZFPDecompressor decompressor;
if (dynField.IsSameType(Handle64()))
{
vtkm::cont::ArrayHandle<Scalar> handle;
const vtkm::Id size = dynField.Cast<Handle64>().GetNumberOfValues();
handle.Allocate(size);
auto fPortal = dynField.Cast<Handle64>().GetPortalControl();
auto hPortal = handle.GetPortalControl();
for (vtkm::Id i = 0; i < size; ++i)
{
hPortal.Set(i, static_cast<Scalar>(fPortal.Get(i)));
}
auto compressed = compressor.Compress(handle, rate, dims);
vtkm::cont::ArrayHandle<Scalar> decoded;
decompressor.Decompress(compressed, decoded, rate, dims);
auto oport = decoded.GetPortalConstControl();
for (int i = 0; i < 4; i++)
{
std::cout << oport.Get(i) << " " << fPortal.Get(i) << " " << oport.Get(i) - fPortal.Get(i)
<< std::endl;
}
}
}
void TestZFP()
{
Test3D<vtkm::Float64>(4);
Test2D<vtkm::Float64>(4);
Test1D<vtkm::Float64>(4);
//Test3D<vtkm::Float32>(4);
//Test3D<vtkm::Int64>(4);
//Test3D<vtkm::Int32>(4);
}
int UnitTestZFPCompressor(int argc, char* argv[])
{
return vtkm::cont::testing::Testing::Run(TestZFP, argc, argv);
}

@ -0,0 +1,40 @@
##============================================================================
## Copyright (c) Kitware, Inc.
## All rights reserved.
## See LICENSE.txt for details.
## This software is distributed WITHOUT ANY WARRANTY; without even
## the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
## PURPOSE. See the above copyright notice for more information.
##
## Copyright 2016 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
## Copyright 2016 UT-Battelle, LLC.
## Copyright 2016 Los Alamos National Security.
##
## Under the terms of Contract DE-NA0003525 with NTESS,
## the U.S. Government retains certain rights in this software.
##
## Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
## Laboratory (LANL), the U.S. Government retains certain rights in
## this software.
##============================================================================
set(headers
ZFPBlockReader.h
ZFPBlockWriter.h
ZFPCodec.h
ZFPDecode.h
ZFPDecode1.h
ZFPDecode2.h
ZFPDecode3.h
ZFPEncode.h
ZFPEncode1.h
ZFPEncode2.h
ZFPEncode3.h
ZFPFunctions.h
ZFPStructs.h
ZFPTools.h
ZFPTypeInfo.h
)
#-----------------------------------------------------------------------------
vtkm_declare_headers(${headers})

@ -0,0 +1,124 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_block_reader_h
#define vtk_m_worklet_zfp_block_reader_h
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
using Word = vtkm::UInt64;
template <vtkm::Int32 block_size, typename WordsPortalType>
struct BlockReader
{
const WordsPortalType& Words;
const vtkm::Int32 m_maxbits;
int m_block_idx;
vtkm::Int32 m_current_bit;
vtkm::Id Index;
Word m_buffer;
const vtkm::Id MaxIndex;
VTKM_EXEC
BlockReader(const WordsPortalType& words, const int& maxbits, const int& block_idx)
: Words(words)
, m_maxbits(maxbits)
, MaxIndex(words.GetNumberOfValues() - 1)
{
Index = vtkm::Id(static_cast<size_t>(block_idx * maxbits) / (sizeof(Word) * 8));
m_buffer = static_cast<Word>(Words.Get(Index));
m_current_bit = (block_idx * maxbits) % vtkm::Int32((sizeof(Word) * 8));
m_buffer >>= m_current_bit;
m_block_idx = block_idx;
//std::cout<<"Reader index "<<Index<<"\n";
//print_bits(m_buffer);
//print_bits(Words.Get(Index));
}
inline VTKM_EXEC unsigned int read_bit()
{
vtkm::UInt32 bit = vtkm::UInt32(m_buffer) & 1u;
++m_current_bit;
m_buffer >>= 1;
// handle moving into next word
if (m_current_bit >= vtkm::Int32(sizeof(Word) * 8))
{
m_current_bit = 0;
++Index;
if (Index > MaxIndex)
return false;
m_buffer = static_cast<Word>(Words.Get(Index));
}
return bit;
}
// note this assumes that n_bits is <= 64
inline VTKM_EXEC vtkm::UInt64 read_bits(const int& n_bits)
{
vtkm::UInt64 bits;
// rem bits will always be positive
vtkm::Int32 rem_bits = vtkm::Int32(sizeof(Word) * 8) - m_current_bit;
vtkm::Int32 first_read = vtkm::Min(rem_bits, n_bits);
// first mask
Word mask = ((Word)1 << ((first_read))) - 1;
bits = m_buffer & mask;
m_buffer >>= n_bits;
m_current_bit += first_read;
vtkm::Int32 next_read = 0;
if (n_bits >= rem_bits)
{
m_current_bit = 0;
// just read in 0s if someone asks for more bits past the end of the array.
// not sure what the best way to deal with this i
Index = vtkm::Min(MaxIndex, Index + 1);
m_buffer = static_cast<Word>(Words.Get(Index));
next_read = n_bits - first_read;
}
// this is basically a no-op when first read constained
// all the bits. TODO: if we have aligned reads, this could
// be a conditional without divergence
mask = ((Word)1 << ((next_read))) - 1;
bits += (m_buffer & mask) << first_read;
m_buffer >>= next_read;
m_current_bit += next_read;
return bits;
}
private:
VTKM_EXEC BlockReader() {}
}; // block reader
} // namespace zfp
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_type_info_h

@ -0,0 +1,188 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_block_writer_h
#define vtk_m_worklet_zfp_block_writer_h
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
using Word = vtkm::UInt64;
template <int block_size, typename AtomicPortalType>
struct BlockWriter
{
union UIntInt {
vtkm::UInt64 uintpart;
vtkm::Int64 intpart;
};
vtkm::Id m_word_index;
vtkm::Int32 m_start_bit;
vtkm::Int32 m_current_bit;
const int m_maxbits;
AtomicPortalType& Portal;
//int debug_index;
VTKM_EXEC BlockWriter(AtomicPortalType& portal, const int& maxbits, const vtkm::Id& block_idx)
: m_current_bit(0)
, m_maxbits(maxbits)
, Portal(portal)
{
m_word_index = (block_idx * maxbits) / vtkm::Int32(sizeof(Word) * 8);
// debug_index = m_word_index;
//std::cout<<"** Block "<<block_idx<<" start "<<m_word_index<<"\n";
m_start_bit = vtkm::Int32((block_idx * maxbits) % vtkm::Int32(sizeof(Word) * 8));
}
template <typename T>
void print_bits(T bits)
{
const int bit_size = sizeof(T) * 8;
for (int i = bit_size - 1; i >= 0; --i)
{
T one = 1;
T mask = one << i;
int val = (bits & mask) >> i;
printf("%d", val);
}
printf("\n");
}
void print()
{
//vtkm::Int64 v = Portal.Add(debug_index,0);
//std::cout<<"current bit "<<m_current_bit<<" debug_index "<<debug_index<<" ";
//print_bits(*reinterpret_cast<vtkm::UInt64*>(&v));
}
// void print(int index)
// {
// vtkm::Int64 v = Portal.Add(index, 0);
// //print_bits(*reinterpret_cast<vtkm::UInt64*>(&v));
// }
inline VTKM_EXEC void Add(const vtkm::Id index, Word& value)
{
UIntInt newval;
UIntInt old;
(void)old;
newval.uintpart = value;
Portal.Add(index, newval.intpart);
//old.uintpart = 0;
//UIntInt expected;
//expected.uintpart = newval.uintpart;
//while(old.uintpart != expected.uintpart)
//{
// expected.uintpart = old.uintpart + newval.uintpart;
// old.intpart = Portal.CompareAndSwap(index, expected.intpart, old.intpart);
//}
}
inline VTKM_EXEC
//void write_bits(const unsigned int &bits, const uint &n_bits, const uint &bit_offset)
vtkm::UInt64
write_bits(const vtkm::UInt64& bits, const unsigned int& n_bits)
{
//std::cout<<"write nbits "<<n_bits<<" "<<m_current_bit<<"\n";
//bool print = m_word_index == 0 && m_start_bit == 0;
const int wbits = sizeof(Word) * 8;
//if(bits == 0) { printf("no\n"); return;}
//uint seg_start = (m_start_bit + bit_offset) % wbits;
//int write_index = m_word_index + (m_start_bit + bit_offset) / wbits;
unsigned int seg_start = (m_start_bit + m_current_bit) % wbits;
vtkm::Id write_index = m_word_index;
write_index += vtkm::Id((m_start_bit + m_current_bit) / wbits);
unsigned int seg_end = seg_start + n_bits - 1;
//int write_index = m_word_index;
unsigned int shift = seg_start;
// we may be asked to write less bits than exist in 'bits'
// so we have to make sure that anything after n is zero.
// If this does not happen, then we may write into a zfp
// block not at the specified index
// uint zero_shift = sizeof(Word) * 8 - n_bits;
Word left = (bits >> n_bits) << n_bits;
Word b = bits - left;
Word add = b << shift;
Add(write_index, add);
//debug_index = write_index;
//if(write_index == 0)
//{
// std::cout<<"*******\n";
// std::cout<<"Current bit "<<m_current_bit<<" writing ";
// print_bits(add);
// print();
// std::cout<<"*******\n";
//}
// n_bits straddles the word boundary
bool straddle = seg_start < sizeof(Word) * 8 && seg_end >= sizeof(Word) * 8;
if (straddle)
{
Word rem = b >> (sizeof(Word) * 8 - shift);
Add(write_index + 1, rem);
//std::cout<<"======\n";
//print_bits(rem);
//std::cout<<"======\n";
// printf("Straddle "); print_bits(rem);
//debug_index = write_index +1;
}
m_current_bit += n_bits;
return bits >> (Word)n_bits;
}
// TODO: optimize
vtkm::UInt32 VTKM_EXEC write_bit(const unsigned int& bit)
{
//bool print = m_word_index == 0 && m_start_bit == 0;
const int wbits = sizeof(Word) * 8;
//if(bits == 0) { printf("no\n"); return;}
//uint seg_start = (m_start_bit + bit_offset) % wbits;
//int write_index = m_word_index + (m_start_bit + bit_offset) / wbits;
unsigned int seg_start = (m_start_bit + m_current_bit) % wbits;
vtkm::Id write_index = m_word_index;
write_index += vtkm::Id((m_start_bit + m_current_bit) / wbits);
//uint seg_end = seg_start;
//int write_index = m_word_index;
unsigned int shift = seg_start;
// we may be asked to write less bits than exist in 'bits'
// so we have to make sure that anything after n is zero.
// If this does not happen, then we may write into a zfp
// block not at the specified index
// uint zero_shift = sizeof(Word) * 8 - n_bits;
Word add = (Word)bit << shift;
Add(write_index, add);
m_current_bit += 1;
return bit;
}
};
} // namespace zfp
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_type_info_h

189
vtkm/worklet/zfp/ZFPCodec.h Normal file

@ -0,0 +1,189 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_codec_h
#define vtk_m_worklet_zfp_codec_h
#define index3(i, j, k) ((i) + 4 * ((j) + 4 * (k)))
#define index2(i, j) ((i) + 4 * (j))
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <vtkm::Int32 BlockSize>
struct ZFPCodec;
template <>
struct ZFPCodec<4>
{
VTKM_EXEC_CONT ZFPCodec() {}
VTKM_EXEC vtkm::UInt8 CodecLookup(vtkm::Int32 x) const
{
VTKM_STATIC_CONSTEXPR_ARRAY vtkm::UInt8 perm_1[4] = { 0, 1, 2, 3 };
return perm_1[x];
}
};
template <>
struct ZFPCodec<16>
{
VTKM_EXEC_CONT ZFPCodec() {}
VTKM_EXEC vtkm::UInt8 CodecLookup(vtkm::Int32 x) const
{
/* order coefficients (i, j) by i + j, then i^2 + j^2 */
VTKM_STATIC_CONSTEXPR_ARRAY vtkm::UInt8 perm_2[16] = {
index2(0, 0), /* 0 : 0 */
index2(1, 0), /* 1 : 1 */
index2(0, 1), /* 2 : 1 */
index2(1, 1), /* 3 : 2 */
index2(2, 0), /* 4 : 2 */
index2(0, 2), /* 5 : 2 */
index2(2, 1), /* 6 : 3 */
index2(1, 2), /* 7 : 3 */
index2(3, 0), /* 8 : 3 */
index2(0, 3), /* 9 : 3 */
index2(2, 2), /* 10 : 4 */
index2(3, 1), /* 11 : 4 */
index2(1, 3), /* 12 : 4 */
index2(3, 2), /* 13 : 5 */
index2(2, 3), /* 14 : 5 */
index2(3, 3), /* 15 : 6 */
};
return perm_2[x];
}
};
template <>
struct ZFPCodec<64>
{
VTKM_EXEC_CONT ZFPCodec() {}
VTKM_EXEC vtkm::UInt8 CodecLookup(vtkm::Int32 x) const
{
/* order coefficients (i, j, k) by i + j + k, then i^2 + j^2 + k^2 */
VTKM_STATIC_CONSTEXPR_ARRAY vtkm::UInt8 perm_3[64] = {
index3(0, 0, 0), /* 0 : 0 */
index3(1, 0, 0), /* 1 : 1 */
index3(0, 1, 0), /* 2 : 1 */
index3(0, 0, 1), /* 3 : 1 */
index3(0, 1, 1), /* 4 : 2 */
index3(1, 0, 1), /* 5 : 2 */
index3(1, 1, 0), /* 6 : 2 */
index3(2, 0, 0), /* 7 : 2 */
index3(0, 2, 0), /* 8 : 2 */
index3(0, 0, 2), /* 9 : 2 */
index3(1, 1, 1), /* 10 : 3 */
index3(2, 1, 0), /* 11 : 3 */
index3(2, 0, 1), /* 12 : 3 */
index3(0, 2, 1), /* 13 : 3 */
index3(1, 2, 0), /* 14 : 3 */
index3(1, 0, 2), /* 15 : 3 */
index3(0, 1, 2), /* 16 : 3 */
index3(3, 0, 0), /* 17 : 3 */
index3(0, 3, 0), /* 18 : 3 */
index3(0, 0, 3), /* 19 : 3 */
index3(2, 1, 1), /* 20 : 4 */
index3(1, 2, 1), /* 21 : 4 */
index3(1, 1, 2), /* 22 : 4 */
index3(0, 2, 2), /* 23 : 4 */
index3(2, 0, 2), /* 24 : 4 */
index3(2, 2, 0), /* 25 : 4 */
index3(3, 1, 0), /* 26 : 4 */
index3(3, 0, 1), /* 27 : 4 */
index3(0, 3, 1), /* 28 : 4 */
index3(1, 3, 0), /* 29 : 4 */
index3(1, 0, 3), /* 30 : 4 */
index3(0, 1, 3), /* 31 : 4 */
index3(1, 2, 2), /* 32 : 5 */
index3(2, 1, 2), /* 33 : 5 */
index3(2, 2, 1), /* 34 : 5 */
index3(3, 1, 1), /* 35 : 5 */
index3(1, 3, 1), /* 36 : 5 */
index3(1, 1, 3), /* 37 : 5 */
index3(3, 2, 0), /* 38 : 5 */
index3(3, 0, 2), /* 39 : 5 */
index3(0, 3, 2), /* 40 : 5 */
index3(2, 3, 0), /* 41 : 5 */
index3(2, 0, 3), /* 42 : 5 */
index3(0, 2, 3), /* 43 : 5 */
index3(2, 2, 2), /* 44 : 6 */
index3(3, 2, 1), /* 45 : 6 */
index3(3, 1, 2), /* 46 : 6 */
index3(1, 3, 2), /* 47 : 6 */
index3(2, 3, 1), /* 48 : 6 */
index3(2, 1, 3), /* 49 : 6 */
index3(1, 2, 3), /* 50 : 6 */
index3(0, 3, 3), /* 51 : 6 */
index3(3, 0, 3), /* 52 : 6 */
index3(3, 3, 0), /* 53 : 6 */
index3(3, 2, 2), /* 54 : 7 */
index3(2, 3, 2), /* 55 : 7 */
index3(2, 2, 3), /* 56 : 7 */
index3(1, 3, 3), /* 57 : 7 */
index3(3, 1, 3), /* 58 : 7 */
index3(3, 3, 1), /* 59 : 7 */
index3(2, 3, 3), /* 60 : 8 */
index3(3, 2, 3), /* 61 : 8 */
index3(3, 3, 2), /* 62 : 8 */
index3(3, 3, 3), /* 63 : 9 */
};
return perm_3[x];
}
};
#undef index3
#undef index2
}
}
} // namespace vtkm::worklet::zfp
#endif

@ -0,0 +1,298 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_decode_h
#define vtk_m_worklet_zfp_decode_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/zfp/ZFPBlockReader.h>
#include <vtkm/worklet/zfp/ZFPCodec.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Int, typename Scalar>
inline VTKM_EXEC Scalar dequantize(const Int& x, const int& e);
template <>
inline VTKM_EXEC vtkm::Float64 dequantize<vtkm::Int64, vtkm::Float64>(const vtkm::Int64& x,
const vtkm::Int32& e)
{
return vtkm::Ldexp((vtkm::Float64)x, e - (CHAR_BIT * scalar_sizeof<vtkm::Float64>() - 2));
}
template <>
inline VTKM_EXEC vtkm::Float32 dequantize<vtkm::Int32, vtkm::Float32>(const vtkm::Int32& x,
const vtkm::Int32& e)
{
return vtkm::Ldexp((vtkm::Float32)x, e - (CHAR_BIT * scalar_sizeof<vtkm::Float32>() - 2));
}
template <>
inline VTKM_EXEC vtkm::Int32 dequantize<vtkm::Int32, vtkm::Int32>(const vtkm::Int32&,
const vtkm::Int32&)
{
return 1;
}
template <>
inline VTKM_EXEC vtkm::Int64 dequantize<vtkm::Int64, vtkm::Int64>(const vtkm::Int64&,
const vtkm::Int32&)
{
return 1;
}
template <class Int, vtkm::UInt32 s>
VTKM_EXEC static void inv_lift(Int* p)
{
Int x, y, z, w;
x = *p;
p += s;
y = *p;
p += s;
z = *p;
p += s;
w = *p;
p += s;
/*
** non-orthogonal transform
** ( 4 6 -4 -1) (x)
** 1/4 * ( 4 2 4 5) (y)
** ( 4 -2 4 -5) (z)
** ( 4 -6 -4 1) (w)
*/
y += w >> 1;
w -= y >> 1;
y += w;
w <<= 1;
w -= y;
z += x;
x <<= 1;
x -= z;
y += z;
z <<= 1;
z -= y;
w += x;
x <<= 1;
x -= w;
p -= s;
*p = w;
p -= s;
*p = z;
p -= s;
*p = y;
p -= s;
*p = x;
}
template <vtkm::Int64 BlockSize>
struct inv_transform;
template <>
struct inv_transform<64>
{
template <typename Int>
VTKM_EXEC void inv_xform(Int* p)
{
unsigned int x, y, z;
/* transform along z */
for (y = 0; y < 4; y++)
for (x = 0; x < 4; x++)
inv_lift<Int, 16>(p + 1 * x + 4 * y);
/* transform along y */
for (x = 0; x < 4; x++)
for (z = 0; z < 4; z++)
inv_lift<Int, 4>(p + 16 * z + 1 * x);
/* transform along x */
for (z = 0; z < 4; z++)
for (y = 0; y < 4; y++)
inv_lift<Int, 1>(p + 4 * y + 16 * z);
}
};
template <>
struct inv_transform<16>
{
template <typename Int>
VTKM_EXEC void inv_xform(Int* p)
{
for (int x = 0; x < 4; ++x)
{
inv_lift<Int, 4>(p + 1 * x);
}
for (int y = 0; y < 4; ++y)
{
inv_lift<Int, 1>(p + 4 * y);
}
}
};
template <>
struct inv_transform<4>
{
template <typename Int>
VTKM_EXEC void inv_xform(Int* p)
{
inv_lift<Int, 1>(p);
}
};
inline VTKM_EXEC vtkm::Int64 uint2int(vtkm::UInt64 x)
{
return static_cast<vtkm::Int64>((x ^ 0xaaaaaaaaaaaaaaaaull) - 0xaaaaaaaaaaaaaaaaull);
}
inline VTKM_EXEC vtkm::Int32 uint2int(vtkm::UInt32 x)
{
return static_cast<vtkm::Int32>((x ^ 0xaaaaaaaau) - 0xaaaaaaaau);
}
// Note: I die a little inside everytime I write this sort of template
template <vtkm::Int32 BlockSize,
typename PortalType,
template <int Size, typename Portal> class ReaderType,
typename UInt>
VTKM_EXEC void decode_ints(ReaderType<BlockSize, PortalType>& reader,
vtkm::Int32& maxbits,
UInt* data,
const vtkm::Int32 intprec)
{
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
data[i] = 0;
}
vtkm::UInt64 x;
const vtkm::UInt32 kmin = 0;
vtkm::Int32 bits = maxbits;
for (vtkm::UInt32 k = static_cast<vtkm::UInt32>(intprec), n = 0; bits && k-- > kmin;)
{
// read bit plane
vtkm::UInt32 m = vtkm::Min(n, vtkm::UInt32(bits));
bits -= m;
x = reader.read_bits(static_cast<vtkm::Int32>(m));
for (; n < BlockSize && bits && (bits--, reader.read_bit()); x += (Word)1 << n++)
for (; n < (BlockSize - 1) && bits && (bits--, !reader.read_bit()); n++)
;
// deposit bit plane
for (int i = 0; x; i++, x >>= 1)
{
data[i] += (UInt)(x & 1u) << k;
}
}
//for (int i = 0; i < BlockSize; i++)
//{
// std::cout<<"Decomp int "<<i<<" = "<<data[i]<<"\n";
//}
}
template <vtkm::Int32 BlockSize, typename Scalar, typename PortalType>
VTKM_EXEC void zfp_decode(Scalar* fblock,
vtkm::Int32 maxbits,
vtkm::UInt32 blockIdx,
PortalType stream)
{
zfp::BlockReader<BlockSize, PortalType> reader(stream, maxbits, vtkm::Int32(blockIdx));
using Int = typename zfp::zfp_traits<Scalar>::Int;
using UInt = typename zfp::zfp_traits<Scalar>::UInt;
vtkm::UInt32 cont = 1;
if (!zfp::is_int<Scalar>())
{
cont = reader.read_bit();
}
if (cont)
{
vtkm::UInt32 ebits = static_cast<vtkm::UInt32>(zfp::get_ebits<Scalar>()) + 1;
vtkm::UInt32 emax;
if (!zfp::is_int<Scalar>())
{
//std::cout<<"ebits "<<ebits<<"\n";
// read in the shared exponent
//vtkm::UInt64 b = reader.read_bits(ebits - 1);
//print_bits(b);
//std::cout<<"b "<<b<<"\n";
//std::cout<<"ebias "<<zfp::get_ebias<Scalar>()<<"\n";
//emax = vtkm::UInt32(b - zfp::get_ebias<Scalar>());
emax = vtkm::UInt32(reader.read_bits(static_cast<vtkm::Int32>(ebits) - 1));
emax -= static_cast<vtkm::UInt32>(zfp::get_ebias<Scalar>());
//std::cout<<"EMAX "<<emax<<"\n";
}
else
{
// no exponent bits
ebits = 0;
}
maxbits -= ebits;
UInt ublock[BlockSize];
decode_ints<BlockSize>(reader, maxbits, ublock, zfp::get_precision<Scalar>());
Int iblock[BlockSize];
const zfp::ZFPCodec<BlockSize> codec;
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
vtkm::UInt8 idx = codec.CodecLookup(i);
iblock[idx] = uint2int(ublock[i]);
}
//for (int i = 0; i < BlockSize; i++)
//{
// std::cout<<"before xform tid "<<i<<" "<<iblock[i]<<"\n";
//}
inv_transform<BlockSize> trans;
trans.inv_xform(iblock);
//for (int i = 0; i < BlockSize; i++)
//{
// std::cout<<"tid "<<i<<" "<<iblock[i]<<"\n";
//}
Scalar inv_w = dequantize<Int, Scalar>(1, static_cast<vtkm::Int32>(emax));
//std::cout<<"dequantize factor "<<inv_w<<"\n";
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
fblock[i] = inv_w * (Scalar)iblock[i];
}
}
}
}
}
} // namespace vtkm::worklet::zfp
#endif

@ -0,0 +1,132 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_decode1_h
#define vtk_m_worklet_zfp_decode1_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/WorkletMapField.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPDecode.h>
#include <vtkm/worklet/zfp/ZFPFunctions.h>
#include <vtkm/worklet/zfp/ZFPStructs.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void ScatterPartial1(const Scalar* q,
PortalType& scalars,
vtkm::Id offset,
vtkm::Int32 nx)
{
vtkm::Id x;
for (x = 0; x < nx; x++, offset++, q++)
{
scalars.Set(offset, *q);
}
}
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void Scatter1(const Scalar* q, PortalType& scalars, vtkm::Id offset)
{
for (vtkm::Id x = 0; x < 4; x++, ++offset)
{
scalars.Set(offset, *q++);
} // x
}
struct Decode1 : public vtkm::worklet::WorkletMapField
{
protected:
vtkm::Id Dims; // field dims
vtkm::Id PaddedDims; // dims padded to a multiple of zfp block size
vtkm::Id ZFPDims; // zfp block dims
vtkm::UInt32 MaxBits; // bits per zfp block
public:
Decode1(const vtkm::Id dims, const vtkm::Id paddedDims, const vtkm::UInt32 maxbits)
: Dims(dims)
, PaddedDims(paddedDims)
, MaxBits(maxbits)
{
ZFPDims = PaddedDims / 4;
}
using ControlSignature = void(FieldIn<>, WholeArrayOut<>, WholeArrayIn<> bitstream);
using ExecutionSignature = void(_1, _2, _3);
template <typename InputScalarPortal, typename BitstreamPortal>
VTKM_EXEC void operator()(const vtkm::Id blockIdx,
InputScalarPortal& scalars,
BitstreamPortal& stream) const
{
using Scalar = typename InputScalarPortal::ValueType;
constexpr vtkm::Int32 BlockSize = 4;
Scalar fblock[BlockSize];
// clear
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
fblock[i] = static_cast<Scalar>(0);
}
zfp::zfp_decode<BlockSize>(
fblock, vtkm::Int32(MaxBits), static_cast<vtkm::UInt32>(blockIdx), stream);
vtkm::Id zfpBlock;
zfpBlock = blockIdx % ZFPDims;
vtkm::Id logicalStart = zfpBlock * vtkm::Id(4);
//std::cout<<"Block ID "<<blockIdx<<"\n";
//std::cout<<"ZFP Block "<<zfpBlock<<"\n";
//std::cout<<"logicalStart Start "<<logicalStart<<"\n";
// get the offset into the field
//vtkm::Id offset = (zfpBlock[2]*4*ZFPDims[1] + zfpBlock[1] * 4)*ZFPDims[0] * 4 + zfpBlock[0] * 4;
//std::cout<<"ZFP block offset "<<offset<<"\n";
bool partial = false;
if (logicalStart + 4 > Dims)
partial = true;
//std::cout<<"Dims "<<Dims<<"\n";
if (partial)
{
const vtkm::Int32 nx =
logicalStart + 4 > Dims ? vtkm::Int32(Dims - logicalStart) : vtkm::Int32(4);
//std::cout<<"Partial block "<<logicalStart<<" offset "<<offset<<"\n";
//std::cout<<"Nx "<<nx<<" "<<ny<<" "<<nz<<"\n";
ScatterPartial1(fblock, scalars, logicalStart, nx);
}
else
{
//std::cout<<"FULL block "<<zfpBlock<<"\n";
Scatter1(fblock, scalars, logicalStart);
}
}
};
}
}
} // namespace vtkm::worklet::zfp
#endif

@ -0,0 +1,149 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_decode2_h
#define vtk_m_worklet_zfp_decode2_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/WorkletMapField.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPDecode.h>
#include <vtkm/worklet/zfp/ZFPFunctions.h>
#include <vtkm/worklet/zfp/ZFPStructs.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void ScatterPartial2(const Scalar* q,
PortalType& scalars,
const vtkm::Id2 dims,
vtkm::Id offset,
vtkm::Int32 nx,
vtkm::Int32 ny)
{
vtkm::Id x, y;
for (y = 0; y < ny; y++, offset += dims[0] - nx, q += 4 - nx)
{
for (x = 0; x < nx; x++, offset++, q++)
{
scalars.Set(offset, *q);
}
}
}
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void Scatter2(const Scalar* q,
PortalType& scalars,
const vtkm::Id2 dims,
vtkm::Id offset)
{
for (vtkm::Id y = 0; y < 4; y++, offset += dims[0] - 4)
{
for (vtkm::Id x = 0; x < 4; x++, ++offset)
{
scalars.Set(offset, *q++);
} // x
} // y
}
struct Decode2 : public vtkm::worklet::WorkletMapField
{
protected:
vtkm::Id2 Dims; // field dims
vtkm::Id2 PaddedDims; // dims padded to a multiple of zfp block size
vtkm::Id2 ZFPDims; // zfp block dims
vtkm::UInt32 MaxBits; // bits per zfp block
public:
Decode2(const vtkm::Id2 dims, const vtkm::Id2 paddedDims, const vtkm::UInt32 maxbits)
: Dims(dims)
, PaddedDims(paddedDims)
, MaxBits(maxbits)
{
ZFPDims[0] = PaddedDims[0] / 4;
ZFPDims[1] = PaddedDims[1] / 4;
}
using ControlSignature = void(FieldIn<>, WholeArrayOut<>, WholeArrayIn<> bitstream);
using ExecutionSignature = void(_1, _2, _3);
template <typename InputScalarPortal, typename BitstreamPortal>
VTKM_EXEC void operator()(const vtkm::Id blockIdx,
InputScalarPortal& scalars,
BitstreamPortal& stream) const
{
using Scalar = typename InputScalarPortal::ValueType;
constexpr vtkm::Int32 BlockSize = 16;
Scalar fblock[BlockSize];
// clear
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
fblock[i] = static_cast<Scalar>(0);
}
zfp::zfp_decode<BlockSize>(
fblock, vtkm::Int32(MaxBits), static_cast<vtkm::UInt32>(blockIdx), stream);
vtkm::Id2 zfpBlock;
zfpBlock[0] = blockIdx % ZFPDims[0];
zfpBlock[1] = (blockIdx / ZFPDims[0]) % ZFPDims[1];
vtkm::Id2 logicalStart = zfpBlock * vtkm::Id(4);
//std::cout<<"Block ID "<<blockIdx<<"\n";
//std::cout<<"ZFP Block "<<zfpBlock<<"\n";
//std::cout<<"logicalStart Start "<<logicalStart<<"\n";
// get the offset into the field
//vtkm::Id offset = (zfpBlock[2]*4*ZFPDims[1] + zfpBlock[1] * 4)*ZFPDims[0] * 4 + zfpBlock[0] * 4;
vtkm::Id offset = logicalStart[0] + logicalStart[1] * Dims[0];
//std::cout<<"ZFP block offset "<<offset<<"\n";
bool partial = false;
if (logicalStart[0] + 4 > Dims[0])
partial = true;
if (logicalStart[1] + 4 > Dims[1])
partial = true;
//std::cout<<"Dims "<<Dims<<"\n";
if (partial)
{
const vtkm::Int32 nx =
logicalStart[0] + 4 > Dims[0] ? vtkm::Int32(Dims[0] - logicalStart[0]) : vtkm::Int32(4);
const vtkm::Int32 ny =
logicalStart[1] + 4 > Dims[1] ? vtkm::Int32(Dims[1] - logicalStart[1]) : vtkm::Int32(4);
//std::cout<<"Partial block "<<logicalStart<<" offset "<<offset<<"\n";
//std::cout<<"Nx "<<nx<<" "<<ny<<" "<<nz<<"\n";
ScatterPartial2(fblock, scalars, Dims, offset, nx, ny);
}
else
{
//std::cout<<"FULL block "<<zfpBlock<<"\n";
Scatter2(fblock, scalars, Dims, offset);
}
}
};
}
}
} // namespace vtkm::worklet::zfp
#endif

@ -0,0 +1,169 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_decode3_h
#define vtk_m_worklet_zfp_decode3_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/WorkletMapField.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPDecode.h>
#include <vtkm/worklet/zfp/ZFPFunctions.h>
#include <vtkm/worklet/zfp/ZFPStructs.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void ScatterPartial3(const Scalar* q,
PortalType& scalars,
const vtkm::Id3 dims,
vtkm::Id offset,
vtkm::Int32 nx,
vtkm::Int32 ny,
vtkm::Int32 nz)
{
vtkm::Id x, y, z;
for (z = 0; z < nz; z++, offset += dims[0] * dims[1] - ny * dims[0], q += 4 * (4 - ny))
{
for (y = 0; y < ny; y++, offset += dims[0] - nx, q += 4 - nx)
{
for (x = 0; x < nx; x++, offset++, q++)
{
scalars.Set(offset, *q);
}
}
}
}
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void Scatter3(const Scalar* q,
PortalType& scalars,
const vtkm::Id3 dims,
vtkm::Id offset)
{
for (vtkm::Id z = 0; z < 4; z++, offset += dims[0] * dims[1] - 4 * dims[0])
{
for (vtkm::Id y = 0; y < 4; y++, offset += dims[0] - 4)
{
for (vtkm::Id x = 0; x < 4; x++, ++offset)
{
scalars.Set(offset, *q++);
} // x
} // y
} // z
}
struct Decode3 : public vtkm::worklet::WorkletMapField
{
protected:
vtkm::Id3 Dims; // field dims
vtkm::Id3 PaddedDims; // dims padded to a multiple of zfp block size
vtkm::Id3 ZFPDims; // zfp block dims
vtkm::UInt32 MaxBits; // bits per zfp block
public:
Decode3(const vtkm::Id3 dims, const vtkm::Id3 paddedDims, const vtkm::UInt32 maxbits)
: Dims(dims)
, PaddedDims(paddedDims)
, MaxBits(maxbits)
{
ZFPDims[0] = PaddedDims[0] / 4;
ZFPDims[1] = PaddedDims[1] / 4;
ZFPDims[2] = PaddedDims[2] / 4;
}
using ControlSignature = void(FieldIn<>, WholeArrayOut<>, WholeArrayIn<> bitstream);
using ExecutionSignature = void(_1, _2, _3);
template <typename InputScalarPortal, typename BitstreamPortal>
VTKM_EXEC void operator()(const vtkm::Id blockIdx,
InputScalarPortal& scalars,
BitstreamPortal& stream) const
{
using Scalar = typename InputScalarPortal::ValueType;
constexpr vtkm::Int32 BlockSize = 64;
Scalar fblock[BlockSize];
// clear
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
fblock[i] = static_cast<Scalar>(0);
}
zfp::zfp_decode<BlockSize>(
fblock, vtkm::Int32(MaxBits), static_cast<vtkm::UInt32>(blockIdx), stream);
//for(int i = 0; i < BlockSize; ++i)
//{
// std::cout<<" "<<fblock[i];
//}
//std::cout<<"\n";
vtkm::Id3 zfpBlock;
zfpBlock[0] = blockIdx % ZFPDims[0];
zfpBlock[1] = (blockIdx / ZFPDims[0]) % ZFPDims[1];
zfpBlock[2] = blockIdx / (ZFPDims[0] * ZFPDims[1]);
vtkm::Id3 logicalStart = zfpBlock * vtkm::Id(4);
//std::cout<<"Block ID "<<blockIdx<<"\n";
//std::cout<<"ZFP Block "<<zfpBlock<<"\n";
//std::cout<<"logicalStart Start "<<logicalStart<<"\n";
// get the offset into the field
//vtkm::Id offset = (zfpBlock[2]*4*ZFPDims[1] + zfpBlock[1] * 4)*ZFPDims[0] * 4 + zfpBlock[0] * 4;
vtkm::Id offset = (logicalStart[2] * Dims[1] + logicalStart[1]) * Dims[0] + logicalStart[0];
//std::cout<<"ZFP block offset "<<offset<<"\n";
bool partial = false;
if (logicalStart[0] + 4 > Dims[0])
partial = true;
if (logicalStart[1] + 4 > Dims[1])
partial = true;
if (logicalStart[2] + 4 > Dims[2])
partial = true;
//std::cout<<"Dims "<<Dims<<"\n";
if (partial)
{
const vtkm::Int32 nx =
logicalStart[0] + 4 > Dims[0] ? vtkm::Int32(Dims[0] - logicalStart[0]) : vtkm::Int32(4);
const vtkm::Int32 ny =
logicalStart[1] + 4 > Dims[1] ? vtkm::Int32(Dims[1] - logicalStart[1]) : vtkm::Int32(4);
const vtkm::Int32 nz =
logicalStart[2] + 4 > Dims[2] ? vtkm::Int32(Dims[2] - logicalStart[2]) : vtkm::Int32(4);
//std::cout<<"Partial block "<<logicalStart<<" offset "<<offset<<"\n";
//std::cout<<"Nx "<<nx<<" "<<ny<<" "<<nz<<"\n";
ScatterPartial3(fblock, scalars, Dims, offset, nx, ny, nz);
}
else
{
//std::cout<<"FULL block "<<zfpBlock<<"\n";
Scatter3(fblock, scalars, Dims, offset);
}
}
};
}
}
} // namespace vtkm::worklet::zfp
#endif

@ -0,0 +1,399 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_encode_h
#define vtk_m_worklet_zfp_encode_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPCodec.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Scalar>
VTKM_EXEC void PadBlock(Scalar* p, vtkm::UInt32 n, vtkm::UInt32 s)
{
switch (n)
{
case 0:
p[0 * s] = 0;
/* FALLTHROUGH */
case 1:
p[1 * s] = p[0 * s];
/* FALLTHROUGH */
case 2:
p[2 * s] = p[1 * s];
/* FALLTHROUGH */
case 3:
p[3 * s] = p[0 * s];
/* FALLTHROUGH */
default:
break;
}
}
template <vtkm::Int32 N, typename FloatType>
inline VTKM_EXEC vtkm::Int32 MaxExponent(const FloatType* vals)
{
FloatType maxVal = 0;
for (vtkm::Int32 i = 0; i < N; ++i)
{
maxVal = vtkm::Max(maxVal, vtkm::Abs(vals[i]));
}
if (maxVal > 0)
{
vtkm::Int32 exponent;
vtkm::Frexp(maxVal, &exponent);
/* clamp exponent in case x is denormal */
return vtkm::Max(exponent, 1 - get_ebias<FloatType>());
}
return -get_ebias<FloatType>();
}
// maximum number of bit planes to encode
inline VTKM_EXEC vtkm::Int32 precision(vtkm::Int32 maxexp, vtkm::Int32 maxprec, vtkm::Int32 minexp)
{
return vtkm::Min(maxprec, vtkm::Max(0, maxexp - minexp + 8));
}
template <typename Scalar>
inline VTKM_EXEC Scalar quantize(Scalar x, vtkm::Int32 e)
{
return vtkm::Ldexp(x, (CHAR_BIT * (vtkm::Int32)sizeof(Scalar) - 2) - e);
}
template <typename Int, typename Scalar, vtkm::Int32 BlockSize>
inline VTKM_EXEC void fwd_cast(Int* iblock, const Scalar* fblock, vtkm::Int32 emax)
{
Scalar s = quantize<Scalar>(1, emax);
//std::cout<<"EMAX "<<emax<<" q "<<s<<"\n";
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
iblock[i] = static_cast<Int>(s * fblock[i]);
//std::cout<<i<<" f = "<<fblock[i]<<" i = "<<(vtkm::UInt64)iblock[i]<<"\n";
}
}
template <typename Int, vtkm::Int32 S>
inline VTKM_EXEC void fwd_lift(Int* p)
{
Int x, y, z, w;
x = *p;
p += S;
y = *p;
p += S;
z = *p;
p += S;
w = *p;
p += S;
/*
** non-orthogonal transform
** ( 4 4 4 4) (x)
** 1/16 * ( 5 1 -1 -5) (y)
** (-4 4 4 -4) (z)
** (-2 6 -6 2) (w)
*/
x += w;
x >>= 1;
w -= x;
z += y;
z >>= 1;
y -= z;
x += z;
x >>= 1;
z -= x;
w += y;
w >>= 1;
y -= w;
w += y >> 1;
y -= w >> 1;
p -= S;
*p = w;
p -= S;
*p = z;
p -= S;
*p = y;
p -= S;
*p = x;
}
template <typename Int, typename UInt>
inline VTKM_EXEC UInt int2uint(const Int x);
template <>
inline VTKM_EXEC vtkm::UInt64 int2uint<vtkm::Int64, vtkm::UInt64>(const vtkm::Int64 x)
{
return (static_cast<vtkm::UInt64>(x) + (vtkm::UInt64)0xaaaaaaaaaaaaaaaaull) ^
(vtkm::UInt64)0xaaaaaaaaaaaaaaaaull;
}
template <>
inline VTKM_EXEC vtkm::UInt32 int2uint<vtkm::Int32, vtkm::UInt32>(const vtkm::Int32 x)
{
return (static_cast<vtkm::UInt32>(x) + (vtkm::UInt32)0xaaaaaaaau) ^ (vtkm::UInt32)0xaaaaaaaau;
}
template <typename UInt, typename Int, vtkm::Int32 BlockSize>
inline VTKM_EXEC void fwd_order(UInt* ublock, const Int* iblock)
{
const zfp::ZFPCodec<BlockSize> codec;
for (vtkm::Int32 i = 0; i < BlockSize; ++i)
{
vtkm::UInt8 idx = codec.CodecLookup(i);
ublock[i] = int2uint<Int, UInt>(iblock[idx]);
}
}
template <typename Int, vtkm::Int32 BlockSize>
inline VTKM_EXEC void fwd_xform(Int* p);
template <>
inline VTKM_EXEC void fwd_xform<vtkm::Int64, 64>(vtkm::Int64* p)
{
vtkm::UInt32 x, y, z;
/* transform along x */
for (z = 0; z < 4; z++)
for (y = 0; y < 4; y++)
fwd_lift<vtkm::Int64, 1>(p + 4 * y + 16 * z);
/* transform along y */
for (x = 0; x < 4; x++)
for (z = 0; z < 4; z++)
fwd_lift<vtkm::Int64, 4>(p + 16 * z + 1 * x);
/* transform along z */
for (y = 0; y < 4; y++)
for (x = 0; x < 4; x++)
fwd_lift<vtkm::Int64, 16>(p + 1 * x + 4 * y);
}
template <>
inline VTKM_EXEC void fwd_xform<vtkm::Int32, 64>(vtkm::Int32* p)
{
vtkm::UInt32 x, y, z;
/* transform along x */
for (z = 0; z < 4; z++)
for (y = 0; y < 4; y++)
fwd_lift<vtkm::Int32, 1>(p + 4 * y + 16 * z);
/* transform along y */
for (x = 0; x < 4; x++)
for (z = 0; z < 4; z++)
fwd_lift<vtkm::Int32, 4>(p + 16 * z + 1 * x);
/* transform along z */
for (y = 0; y < 4; y++)
for (x = 0; x < 4; x++)
fwd_lift<vtkm::Int32, 16>(p + 1 * x + 4 * y);
}
template <>
inline VTKM_EXEC void fwd_xform<vtkm::Int64, 16>(vtkm::Int64* p)
{
vtkm::UInt32 x, y;
/* transform along x */
for (y = 0; y < 4; y++)
fwd_lift<vtkm::Int64, 1>(p + 4 * y);
/* transform along y */
for (x = 0; x < 4; x++)
fwd_lift<vtkm::Int64, 4>(p + 1 * x);
}
template <>
inline VTKM_EXEC void fwd_xform<vtkm::Int32, 16>(vtkm::Int32* p)
{
vtkm::UInt32 x, y;
/* transform along x */
for (y = 0; y < 4; y++)
fwd_lift<vtkm::Int32, 1>(p + 4 * y);
/* transform along y */
for (x = 0; x < 4; x++)
fwd_lift<vtkm::Int32, 4>(p + 1 * x);
}
template <>
inline VTKM_EXEC void fwd_xform<vtkm::Int64, 4>(vtkm::Int64* p)
{
/* transform along x */
fwd_lift<vtkm::Int64, 1>(p);
}
template <>
inline VTKM_EXEC void fwd_xform<vtkm::Int32, 4>(vtkm::Int32* p)
{
/* transform along x */
fwd_lift<vtkm::Int32, 1>(p);
}
template <vtkm::Int32 BlockSize, typename PortalType, typename Int>
VTKM_EXEC void encode_block(BlockWriter<BlockSize, PortalType>& stream,
vtkm::Int32 maxbits,
vtkm::Int32 maxprec,
Int* iblock)
{
using UInt = typename zfp_traits<Int>::UInt;
fwd_xform<Int, BlockSize>(iblock);
UInt ublock[BlockSize];
fwd_order<UInt, Int, BlockSize>(ublock, iblock);
//for(int i = 0; i < BlockSize; ++i)
//{
// std::cout<<"tid "<<i<<" --> nb "<<ublock[i]<<"\n";
//}
vtkm::UInt32 intprec = CHAR_BIT * (vtkm::UInt32)sizeof(UInt);
vtkm::UInt32 kmin =
intprec > (vtkm::UInt32)maxprec ? intprec - static_cast<vtkm::UInt32>(maxprec) : 0;
vtkm::UInt32 bits = static_cast<vtkm::UInt32>(maxbits);
vtkm::UInt32 i, m;
vtkm::UInt32 n = 0;
vtkm::UInt64 x;
//std::cout<<"Kmin "<<kmin<<"\n";
/* encode one bit plane at a time from MSB to LSB */
for (vtkm::UInt32 k = intprec; bits && k-- > kmin;)
{
/* step 1: extract bit plane #k to x */
x = 0;
for (i = 0; i < BlockSize; i++)
{
x += (vtkm::UInt64)((ublock[i] >> k) & 1u) << i;
}
//std::cout<<"Bit plane "<<x<<"\n";
/* step 2: encode first n bits of bit plane */
m = vtkm::Min(n, bits);
bits -= m;
//std::cout<<"Bits left "<<bits<<" m "<<m<<"\n";
x = stream.write_bits(x, m);
//std::cout<<"Wrote m "<<m<<" bits\n";
//vtkm::UInt32 temp = bits;
//std::cout<<"rem bitplane "<<x<<"\n";
/* step 3: unary run-length encode remainder of bit plane */
for (; n < BlockSize && bits && (bits--, stream.write_bit(!!x)); x >>= 1, n++)
{
//std::cout<<"outer n "<<n<<" bits "<<bits<<"\n";
for (; n < BlockSize - 1 && bits && (bits--, !stream.write_bit(x & 1u)); x >>= 1, n++)
{
//std::cout<<"n "<<n<<" bits "<<bits<<"\n";
}
}
//temp = temp - bits;
//std::cout<<"rem bits "<<bits<<" intprec "<<intprec<<" k "<<k<<" encoded_bits "<<temp<<"\n";
//stream.print();
}
}
template <vtkm::Int32 BlockSize, typename Scalar, typename PortalType>
inline VTKM_EXEC void zfp_encodef(Scalar* fblock,
vtkm::Int32 maxbits,
vtkm::UInt32 blockIdx,
PortalType& stream)
{
using Int = typename zfp::zfp_traits<Scalar>::Int;
zfp::BlockWriter<BlockSize, PortalType> blockWriter(stream, maxbits, vtkm::Id(blockIdx));
vtkm::Int32 emax = zfp::MaxExponent<BlockSize, Scalar>(fblock);
// std::cout<<"EMAX "<<emax<<"\n";
vtkm::Int32 maxprec =
zfp::precision(emax, zfp::get_precision<Scalar>(), zfp::get_min_exp<Scalar>());
vtkm::UInt32 e = vtkm::UInt32(maxprec ? emax + zfp::get_ebias<Scalar>() : 0);
/* encode block only if biased exponent is nonzero */
if (e)
{
const vtkm::UInt32 ebits = vtkm::UInt32(zfp::get_ebits<Scalar>()) + 1;
blockWriter.write_bits(2 * e + 1, ebits);
Int iblock[BlockSize];
zfp::fwd_cast<Int, Scalar, BlockSize>(iblock, fblock, emax);
encode_block<BlockSize>(blockWriter, maxbits - vtkm::Int32(ebits), maxprec, iblock);
}
}
// helpers so we can do partial template instantiation since
// the portal type could be on any backend
template <vtkm::Int32 BlockSize, typename Scalar, typename PortalType>
struct ZFPBlockEncoder
{
};
template <vtkm::Int32 BlockSize, typename PortalType>
struct ZFPBlockEncoder<BlockSize, vtkm::Float32, PortalType>
{
VTKM_EXEC void encode(vtkm::Float32* fblock,
vtkm::Int32 maxbits,
vtkm::UInt32 blockIdx,
PortalType& stream)
{
zfp_encodef<BlockSize>(fblock, maxbits, blockIdx, stream);
}
};
template <vtkm::Int32 BlockSize, typename PortalType>
struct ZFPBlockEncoder<BlockSize, vtkm::Float64, PortalType>
{
VTKM_EXEC void encode(vtkm::Float64* fblock,
vtkm::Int32 maxbits,
vtkm::UInt32 blockIdx,
PortalType& stream)
{
zfp_encodef<BlockSize>(fblock, maxbits, blockIdx, stream);
}
};
template <vtkm::Int32 BlockSize, typename PortalType>
struct ZFPBlockEncoder<BlockSize, vtkm::Int32, PortalType>
{
VTKM_EXEC void encode(vtkm::Int32* fblock,
vtkm::Int32 maxbits,
vtkm::UInt32 blockIdx,
PortalType& stream)
{
using Int = typename zfp::zfp_traits<vtkm::Int32>::Int;
zfp::BlockWriter<BlockSize, PortalType> blockWriter(stream, maxbits, vtkm::Id(blockIdx));
encode_block<BlockSize>(blockWriter, maxbits, get_precision<vtkm::Int32>(), (Int*)fblock);
}
};
template <vtkm::Int32 BlockSize, typename PortalType>
struct ZFPBlockEncoder<BlockSize, vtkm::Int64, PortalType>
{
VTKM_EXEC void encode(vtkm::Int64* fblock,
vtkm::Int32 maxbits,
vtkm::UInt32 blockIdx,
PortalType& stream)
{
using Int = typename zfp::zfp_traits<vtkm::Int64>::Int;
zfp::BlockWriter<BlockSize, PortalType> blockWriter(stream, maxbits, vtkm::Id(blockIdx));
encode_block<BlockSize>(blockWriter, maxbits, get_precision<vtkm::Int64>(), (Int*)fblock);
}
};
}
}
} // namespace vtkm::worklet::zfp
#endif

@ -0,0 +1,155 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_encode1_h
#define vtk_m_worklet_zfp_encode1_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/WorkletMapField.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPEncode.h>
#include <vtkm/worklet/zfp/ZFPFunctions.h>
#include <vtkm/worklet/zfp/ZFPStructs.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void GatherPartial1(Scalar* q,
const PortalType& scalars,
vtkm::Id offset,
int nx,
int sx)
{
vtkm::Id x;
for (x = 0; x < nx; x++, offset += sx)
q[x] = scalars.Get(offset);
PadBlock(q, vtkm::UInt32(nx), 1);
}
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void Gather1(Scalar* fblock, const PortalType& scalars, vtkm::Id offset, int sx)
{
vtkm::Id counter = 0;
for (vtkm::Id x = 0; x < 4; x++, offset += sx)
{
fblock[counter] = scalars.Get(offset);
counter++;
}
}
struct Encode1 : public vtkm::worklet::WorkletMapField
{
protected:
vtkm::Id Dims; // field dims
vtkm::Id PaddedDims; // dims padded to a multiple of zfp block size
vtkm::Id ZFPDims; // zfp block dims
vtkm::UInt32 MaxBits; // bits per zfp block
public:
Encode1(const vtkm::Id dims, const vtkm::Id paddedDims, const vtkm::UInt32 maxbits)
: Dims(dims)
, PaddedDims(paddedDims)
, MaxBits(maxbits)
{
ZFPDims = PaddedDims / 4;
}
using ControlSignature = void(FieldIn<>, WholeArrayIn<>, AtomicArrayInOut<> bitstream);
using ExecutionSignature = void(_1, _2, _3);
template <class InputScalarPortal, typename BitstreamPortal>
VTKM_EXEC void operator()(const vtkm::Id blockIdx,
const InputScalarPortal& scalars,
BitstreamPortal& stream) const
{
using Scalar = typename InputScalarPortal::ValueType;
// typedef unsigned long long int ull;
// typedef long long int ll;
// const ull blockId = blockIdx.x +
// blockIdx.y * gridDim.x +
// gridDim.x * gridDim.y * blockIdx.z;
// // each thread gets a block so the block index is
// // the global thread index
// const uint block_idx = blockId * blockDim.x + threadIdx.x;
// if(block_idx >= tot_blocks)
// {
// // we can't launch the exact number of blocks
// // so just exit if this isn't real
// return;
// }
// uint2 block_dims;
// block_dims.x = padded_dims.x >> 2;
// block_dims.y = padded_dims.y >> 2;
// // logical pos in 3d array
// uint2 block;
// block.x = (block_idx % block_dims.x) * 4;
// block.y = ((block_idx/ block_dims.x) % block_dims.y) * 4;
// const ll offset = (ll)block.x * stride.x + (ll)block.y * stride.y;
vtkm::Id zfpBlock;
zfpBlock = blockIdx % ZFPDims;
vtkm::Id logicalStart = zfpBlock * vtkm::Id(4);
constexpr vtkm::Int32 BlockSize = 4;
Scalar fblock[BlockSize];
// bool partial = false;
// if(block.x + 4 > dims.x) partial = true;
// if(block.y + 4 > dims.y) partial = true;
bool partial = false;
if (logicalStart + 4 > Dims)
partial = true;
if (partial)
{
const vtkm::Int32 nx =
logicalStart + 4 > Dims ? vtkm::Int32(Dims - logicalStart) : vtkm::Int32(4);
GatherPartial1(fblock, scalars, logicalStart, nx, 1);
}
else
{
Gather1(fblock, scalars, logicalStart, 1);
}
//zfp_encode_block<Scalar, ZFP_2D_BLOCK_SIZE>(fblock, maxbits, block_idx, stream);
zfp::ZFPBlockEncoder<BlockSize, Scalar, BitstreamPortal> encoder;
encoder.encode(fblock, static_cast<vtkm::Int32>(MaxBits), vtkm::UInt32(blockIdx), stream);
}
};
}
}
}
#endif

@ -0,0 +1,173 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_encode2_h
#define vtk_m_worklet_zfp_encode2_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/WorkletMapField.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPEncode.h>
#include <vtkm/worklet/zfp/ZFPFunctions.h>
#include <vtkm/worklet/zfp/ZFPStructs.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void GatherPartial2(Scalar* q,
const PortalType& scalars,
vtkm::Id offset,
vtkm::Int32 nx,
vtkm::Int32 ny,
vtkm::Int32 sx,
vtkm::Int32 sy)
{
vtkm::Id x, y;
for (y = 0; y < ny; y++, offset += sy - nx * sx)
{
for (x = 0; x < nx; x++, offset += 1)
q[4 * y + x] = scalars.Get(offset);
PadBlock(q + 4 * y, vtkm::UInt32(nx), 1);
}
for (x = 0; x < 4; x++)
PadBlock(q + x, vtkm::UInt32(ny), 4);
}
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void Gather2(Scalar* fblock,
const PortalType& scalars,
vtkm::Id offset,
int sx,
int sy)
{
vtkm::Id counter = 0;
for (vtkm::Id y = 0; y < 4; y++, offset += sy - 4 * sx)
for (vtkm::Id x = 0; x < 4; x++, offset += sx)
{
fblock[counter] = scalars.Get(offset);
counter++;
}
}
struct Encode2 : public vtkm::worklet::WorkletMapField
{
protected:
vtkm::Id2 Dims; // field dims
vtkm::Id2 PaddedDims; // dims padded to a multiple of zfp block size
vtkm::Id2 ZFPDims; // zfp block dims
vtkm::UInt32 MaxBits; // bits per zfp block
public:
Encode2(const vtkm::Id2 dims, const vtkm::Id2 paddedDims, const vtkm::UInt32 maxbits)
: Dims(dims)
, PaddedDims(paddedDims)
, MaxBits(maxbits)
{
ZFPDims[0] = PaddedDims[0] / 4;
ZFPDims[1] = PaddedDims[1] / 4;
}
using ControlSignature = void(FieldIn<>, WholeArrayIn<>, AtomicArrayInOut<> bitstream);
using ExecutionSignature = void(_1, _2, _3);
template <class InputScalarPortal, typename BitstreamPortal>
VTKM_EXEC void operator()(const vtkm::Id blockIdx,
const InputScalarPortal& scalars,
BitstreamPortal& stream) const
{
using Scalar = typename InputScalarPortal::ValueType;
// typedef unsigned long long int ull;
// typedef long long int ll;
// const ull blockId = blockIdx.x +
// blockIdx.y * gridDim.x +
// gridDim.x * gridDim.y * blockIdx.z;
// // each thread gets a block so the block index is
// // the global thread index
// const uint block_idx = blockId * blockDim.x + threadIdx.x;
// if(block_idx >= tot_blocks)
// {
// // we can't launch the exact number of blocks
// // so just exit if this isn't real
// return;
// }
// uint2 block_dims;
// block_dims.x = padded_dims.x >> 2;
// block_dims.y = padded_dims.y >> 2;
// // logical pos in 3d array
// uint2 block;
// block.x = (block_idx % block_dims.x) * 4;
// block.y = ((block_idx/ block_dims.x) % block_dims.y) * 4;
// const ll offset = (ll)block.x * stride.x + (ll)block.y * stride.y;
vtkm::Id2 zfpBlock;
zfpBlock[0] = blockIdx % ZFPDims[0];
zfpBlock[1] = (blockIdx / ZFPDims[0]) % ZFPDims[1];
vtkm::Id2 logicalStart = zfpBlock * vtkm::Id(4);
vtkm::Id offset = logicalStart[1] * Dims[0] + logicalStart[0];
constexpr vtkm::Int32 BlockSize = 16;
Scalar fblock[BlockSize];
// bool partial = false;
// if(block.x + 4 > dims.x) partial = true;
// if(block.y + 4 > dims.y) partial = true;
bool partial = false;
if (logicalStart[0] + 4 > Dims[0])
partial = true;
if (logicalStart[1] + 4 > Dims[1])
partial = true;
if (partial)
{
const vtkm::Int32 nx =
logicalStart[0] + 4 > Dims[0] ? vtkm::Int32(Dims[0] - logicalStart[0]) : vtkm::Int32(4);
const vtkm::Int32 ny =
logicalStart[1] + 4 > Dims[1] ? vtkm::Int32(Dims[1] - logicalStart[1]) : vtkm::Int32(4);
GatherPartial2(fblock, scalars, offset, nx, ny, 1, static_cast<vtkm::Int32>(Dims[0]));
}
else
{
Gather2(fblock, scalars, offset, 1, static_cast<vtkm::Int32>(Dims[0]));
}
//zfp_encode_block<Scalar, ZFP_2D_BLOCK_SIZE>(fblock, maxbits, block_idx, stream);
zfp::ZFPBlockEncoder<BlockSize, Scalar, BitstreamPortal> encoder;
encoder.encode(fblock, vtkm::Int32(MaxBits), vtkm::UInt32(blockIdx), stream);
}
};
}
}
}
#endif

@ -0,0 +1,182 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_encode3_h
#define vtk_m_worklet_zfp_encode3_h
#include <vtkm/Types.h>
#include <vtkm/internal/ExportMacros.h>
#include <vtkm/worklet/WorkletMapField.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPEncode.h>
#include <vtkm/worklet/zfp/ZFPFunctions.h>
#include <vtkm/worklet/zfp/ZFPStructs.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void GatherPartial3(Scalar* q,
const PortalType& scalars,
const vtkm::Id3 dims,
vtkm::Id offset,
vtkm::Int32 nx,
vtkm::Int32 ny,
vtkm::Int32 nz)
{
vtkm::Id x, y, z;
for (z = 0; z < nz; z++, offset += dims[0] * dims[1] - ny * dims[0])
{
for (y = 0; y < ny; y++, offset += dims[0] - nx)
{
for (x = 0; x < nx; x++, offset += 1)
{
q[16 * z + 4 * y + x] = scalars.Get(offset);
}
PadBlock(q + 16 * z + 4 * y, static_cast<vtkm::UInt32>(nx), 1);
}
for (x = 0; x < 4; x++)
{
PadBlock(q + 16 * z + x, vtkm::UInt32(ny), 4);
}
}
for (y = 0; y < 4; y++)
{
for (x = 0; x < 4; x++)
{
PadBlock(q + 4 * y + x, vtkm::UInt32(nz), 16);
}
}
}
template <typename Scalar, typename PortalType>
VTKM_EXEC inline void Gather3(Scalar* fblock,
const PortalType& scalars,
const vtkm::Id3 dims,
vtkm::Id offset)
{
// TODO: gather partial
vtkm::Id counter = 0;
for (vtkm::Id z = 0; z < 4; z++, offset += dims[0] * dims[1] - 4 * dims[0])
{
for (vtkm::Id y = 0; y < 4; y++, offset += dims[0] - 4)
{
for (vtkm::Id x = 0; x < 4; x++, ++offset)
{
fblock[counter] = scalars.Get(offset);
counter++;
} // x
} // y
} // z
}
struct Encode3 : public vtkm::worklet::WorkletMapField
{
protected:
vtkm::Id3 Dims; // field dims
vtkm::Id3 PaddedDims; // dims padded to a multiple of zfp block size
vtkm::Id3 ZFPDims; // zfp block dims
vtkm::UInt32 MaxBits; // bits per zfp block
public:
Encode3(const vtkm::Id3 dims, const vtkm::Id3 paddedDims, const vtkm::UInt32 maxbits)
: Dims(dims)
, PaddedDims(paddedDims)
, MaxBits(maxbits)
{
ZFPDims[0] = PaddedDims[0] / 4;
ZFPDims[1] = PaddedDims[1] / 4;
ZFPDims[2] = PaddedDims[2] / 4;
}
using ControlSignature = void(FieldIn<>, WholeArrayIn<>, AtomicArrayInOut<> bitstream);
using ExecutionSignature = void(_1, _2, _3);
template <typename InputScalarPortal, typename BitstreamPortal>
VTKM_EXEC void operator()(const vtkm::Id blockIdx,
const InputScalarPortal& scalars,
BitstreamPortal& stream) const
{
using Scalar = typename InputScalarPortal::ValueType;
constexpr vtkm::Int32 BlockSize = 64;
Scalar fblock[BlockSize];
vtkm::Id3 zfpBlock;
zfpBlock[0] = blockIdx % ZFPDims[0];
zfpBlock[1] = (blockIdx / ZFPDims[0]) % ZFPDims[1];
zfpBlock[2] = blockIdx / (ZFPDims[0] * ZFPDims[1]);
vtkm::Id3 logicalStart = zfpBlock * vtkm::Id(4);
//std::cout<<"Block ID "<<blockIdx<<"\n";
//std::cout<<"ZFP Block "<<zfpBlock<<"\n";
//std::cout<<"logicalStart Start "<<logicalStart<<"\n";
// get the offset into the field
//vtkm::Id offset = (zfpBlock[2]*4*ZFPDims[1] + zfpBlock[1] * 4)*ZFPDims[0] * 4 + zfpBlock[0] * 4;
vtkm::Id offset = (logicalStart[2] * Dims[1] + logicalStart[1]) * Dims[0] + logicalStart[0];
//std::cout<<"ZFP block offset "<<offset<<"\n";
bool partial = false;
if (logicalStart[0] + 4 > Dims[0])
partial = true;
if (logicalStart[1] + 4 > Dims[1])
partial = true;
if (logicalStart[2] + 4 > Dims[2])
partial = true;
//std::cout<<"Dims "<<Dims<<"\n";
if (partial)
{
const vtkm::Int32 nx =
logicalStart[0] + 4 > Dims[0] ? vtkm::Int32(Dims[0] - logicalStart[0]) : vtkm::Int32(4);
const vtkm::Int32 ny =
logicalStart[1] + 4 > Dims[1] ? vtkm::Int32(Dims[1] - logicalStart[1]) : vtkm::Int32(4);
const vtkm::Int32 nz =
logicalStart[2] + 4 > Dims[2] ? vtkm::Int32(Dims[2] - logicalStart[2]) : vtkm::Int32(4);
//std::cout<<"Partial block "<<logicalStart<<" offset "<<offset<<"\n";
//std::cout<<"Nx "<<nx<<" "<<ny<<" "<<nz<<"\n";
GatherPartial3(fblock, scalars, Dims, offset, nx, ny, nz);
}
else
{
//std::cout<<"FULL block "<<zfpBlock<<"\n";
Gather3(fblock, scalars, Dims, offset);
}
// for(int i = 0; i < 64; ++i)
// {
// std::cout<<" "<<fblock[i];
// }
// std::cout<<"\n";
zfp::ZFPBlockEncoder<BlockSize, Scalar, BitstreamPortal> encoder;
encoder.encode(fblock, vtkm::Int32(MaxBits), vtkm::UInt32(blockIdx), stream);
//blockWriter.print(0);
//blockWriter.print(1);
}
};
}
}
} // namespace vtkm::worklet::zfp
#endif

@ -0,0 +1,74 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_functions_h
#define vtk_m_worklet_zfp_functions_h
#include <vtkm/Math.h>
#include <vtkm/worklet/zfp/ZFPBlockWriter.h>
#include <vtkm/worklet/zfp/ZFPCodec.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename T>
void print_bits(T bits)
{
const int bit_size = sizeof(T) * 8;
for (int i = bit_size - 1; i >= 0; --i)
{
T one = 1;
T mask = one << i;
int val = (bits & mask) >> T(i);
printf("%d", val);
}
printf("\n");
}
template <typename T>
inline vtkm::UInt32 MinBits(const vtkm::UInt32 bits)
{
return bits;
}
template <>
inline vtkm::UInt32 MinBits<vtkm::Float32>(const vtkm::UInt32 bits)
{
return vtkm::Max(bits, 1 + 8u);
}
template <>
inline vtkm::UInt32 MinBits<vtkm::Float64>(const vtkm::UInt32 bits)
{
return vtkm::Max(bits, 1 + 11u);
}
} // namespace zfp
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_type_info_h

@ -0,0 +1,66 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_structs_h
#define vtk_m_worklet_zfp_structs_h
#define ZFP_MIN_BITS 0 /* minimum number of bits per block */
#define ZFP_MAX_BITS 4171 /* maximum number of bits per block */
#define ZFP_MAX_PREC 64 /* maximum precision supported */
#define ZFP_MIN_EXP -1074 /* minimum floating-point base-2 exponent */
#include <vtkm/worklet/zfp/ZFPFunctions.h>
#include <vtkm/worklet/zfp/ZFPTypeInfo.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
struct ZFPStream
{
vtkm::UInt32 minbits;
vtkm::UInt32 maxbits;
vtkm::UInt32 maxprec;
vtkm::Int32 minexp;
template <typename T>
vtkm::Float64 SetRate(const vtkm::Float64 rate, const vtkm::Int32 dims, T vtkmNotUsed(valueType))
{
vtkm::UInt32 n = 1u << (2 * dims);
vtkm::UInt32 bits = (unsigned int)floor(n * rate + 0.5);
bits = zfp::MinBits<T>(bits);
//if (wra) {
// /* for write random access, round up to next multiple of stream word size */
// bits += (uint)stream_word_bits - 1;
// bits &= ~(stream_word_bits - 1);
//}
minbits = bits;
maxbits = bits;
maxprec = ZFP_MAX_PREC;
minexp = ZFP_MIN_EXP;
return (double)bits / n;
}
};
}
}
} // namespace vtkm::worklet::zfp
#endif

128
vtkm/worklet/zfp/ZFPTools.h Normal file

@ -0,0 +1,128 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_tool_h
#define vtk_m_worklet_zfp_tool_h
#include <vtkm/Math.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/ArrayHandleConstant.h>
#include <vtkm/cont/ArrayHandleCounting.h>
#include <vtkm/cont/AtomicArray.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/zfp/ZFPEncode3.h>
using ZFPWord = vtkm::UInt64;
#include <stdio.h>
namespace vtkm
{
namespace worklet
{
namespace detail
{
class MemTransfer : public vtkm::worklet::WorkletMapField
{
public:
VTKM_CONT
MemTransfer() {}
using ControlSignature = void(FieldIn<>, WholeArrayInOut<>);
using ExecutionSignature = void(_1, _2);
template <typename PortalType>
VTKM_EXEC void operator()(const vtkm::Id id, PortalType& outValue) const
{
(void)id;
(void)outValue;
}
}; //class MemTransfer
inline size_t CalcMem3d(const vtkm::Id3 dims, const vtkm::UInt32 bits_per_block)
{
const size_t vals_per_block = 64;
const size_t size = static_cast<size_t>(dims[0] * dims[1] * dims[2]);
size_t total_blocks = size / vals_per_block;
const size_t bits_per_word = sizeof(ZFPWord) * 8;
const size_t total_bits = bits_per_block * total_blocks;
const size_t alloc_size = total_bits / bits_per_word;
return alloc_size * sizeof(ZFPWord);
}
inline size_t CalcMem2d(const vtkm::Id2 dims, const vtkm::UInt32 bits_per_block)
{
constexpr size_t vals_per_block = 16;
const size_t size = static_cast<size_t>(dims[0] * dims[1]);
size_t total_blocks = size / vals_per_block;
constexpr size_t bits_per_word = sizeof(ZFPWord) * 8;
const size_t total_bits = bits_per_block * total_blocks;
const size_t alloc_size = total_bits / bits_per_word;
return alloc_size * sizeof(ZFPWord);
}
inline size_t CalcMem1d(const vtkm::Id dims, const vtkm::UInt32 bits_per_block)
{
constexpr size_t vals_per_block = 4;
const size_t size = static_cast<size_t>(dims);
size_t total_blocks = size / vals_per_block;
constexpr size_t bits_per_word = sizeof(ZFPWord) * 8;
const size_t total_bits = bits_per_block * total_blocks;
const size_t alloc_size = total_bits / bits_per_word;
return alloc_size * sizeof(ZFPWord);
}
} // namespace detail
template <typename T>
T* GetVTKMPointer(vtkm::cont::ArrayHandle<T>& handle)
{
typedef typename vtkm::cont::ArrayHandle<T> HandleType;
typedef typename HandleType::template ExecutionTypes<vtkm::cont::DeviceAdapterTagSerial>::Portal
PortalType;
typedef typename vtkm::cont::ArrayPortalToIterators<PortalType>::IteratorType IteratorType;
IteratorType iter =
vtkm::cont::ArrayPortalToIterators<PortalType>(handle.GetPortalControl()).GetBegin();
return &(*iter);
}
template <typename T>
void DataDump(vtkm::cont::ArrayHandle<T> handle, std::string fileName)
{
T* ptr = GetVTKMPointer(handle);
vtkm::Id osize = handle.GetNumberOfValues();
FILE* fp = fopen(fileName.c_str(), "wb");
;
if (fp != NULL)
{
fwrite(ptr, sizeof(T), static_cast<size_t>(osize), fp);
}
fclose(fp);
}
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_tools_h

@ -0,0 +1,238 @@
//============================================================================
// Copyright (c) Kitware, Inc.
// All rights reserved.
// See LICENSE.txt for details.
// This software is distributed WITHOUT ANY WARRANTY; without even
// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
// PURPOSE. See the above copyright notice for more information.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef vtk_m_worklet_zfp_type_info_compressor_h
#define vtk_m_worklet_zfp_type_info_compressor_h
#include <vtkm/Math.h>
namespace vtkm
{
namespace worklet
{
namespace zfp
{
template <typename T>
inline VTKM_EXEC int get_ebias();
template <>
inline VTKM_EXEC int get_ebias<double>()
{
return 1023;
}
template <>
inline VTKM_EXEC int get_ebias<float>()
{
return 127;
}
template <>
inline VTKM_EXEC int get_ebias<long long int>()
{
return 0;
}
template <>
inline VTKM_EXEC int get_ebias<int>()
{
return 0;
}
template <typename T>
inline VTKM_EXEC int get_ebits();
template <>
inline VTKM_EXEC int get_ebits<double>()
{
return 11;
}
template <>
inline VTKM_EXEC int get_ebits<float>()
{
return 8;
}
template <>
inline VTKM_EXEC int get_ebits<int>()
{
return 0;
}
template <>
inline VTKM_EXEC int get_ebits<long long int>()
{
return 0;
}
template <typename T>
inline VTKM_EXEC int get_precision();
template <>
inline VTKM_EXEC int get_precision<double>()
{
return 64;
}
template <>
inline VTKM_EXEC int get_precision<long long int>()
{
return 64;
}
template <>
inline VTKM_EXEC int get_precision<float>()
{
return 32;
}
template <>
inline VTKM_EXEC int get_precision<int>()
{
return 32;
}
template <typename T>
inline VTKM_EXEC int get_min_exp();
template <>
inline VTKM_EXEC int get_min_exp<double>()
{
return -1074;
}
template <>
inline VTKM_EXEC int get_min_exp<float>()
{
return -1074;
}
template <>
inline VTKM_EXEC int get_min_exp<long long int>()
{
return 0;
}
template <>
inline VTKM_EXEC int get_min_exp<int>()
{
return 0;
}
template <typename T>
inline VTKM_EXEC int scalar_sizeof();
template <>
inline VTKM_EXEC int scalar_sizeof<double>()
{
return 8;
}
template <>
inline VTKM_EXEC int scalar_sizeof<long long int>()
{
return 8;
}
template <>
inline VTKM_EXEC int scalar_sizeof<float>()
{
return 4;
}
template <>
inline VTKM_EXEC int scalar_sizeof<int>()
{
return 4;
}
template <typename T>
inline VTKM_EXEC bool is_float();
template <>
inline VTKM_EXEC bool is_float<double>()
{
return true;
}
template <>
inline VTKM_EXEC bool is_float<long long int>()
{
return true;
}
template <>
inline VTKM_EXEC bool is_float<float>()
{
return false;
}
template <>
inline VTKM_EXEC bool is_float<int>()
{
return false;
}
template <typename T>
struct zfp_traits;
template <>
struct zfp_traits<double>
{
typedef unsigned long long int UInt;
typedef long long int Int;
};
template <>
struct zfp_traits<long long int>
{
typedef unsigned long long int UInt;
typedef long long int Int;
};
template <>
struct zfp_traits<float>
{
typedef unsigned int UInt;
typedef int Int;
};
template <>
struct zfp_traits<int>
{
typedef unsigned int UInt;
typedef int Int;
};
template <typename T>
inline VTKM_EXEC bool is_int()
{
return false;
}
template <>
inline VTKM_EXEC bool is_int<int>()
{
return true;
}
template <>
inline VTKM_EXEC bool is_int<long long int>()
{
return true;
}
template <int T>
struct block_traits;
template <>
struct block_traits<1>
{
typedef unsigned char PlaneType;
};
template <>
struct block_traits<2>
{
typedef unsigned short PlaneType;
};
} // namespace zfp
} // namespace worklet
} // namespace vtkm
#endif // vtk_m_worklet_zfp_type_info_h