2020-06-17 21:48:09 +00:00
|
|
|
//============================================================================
|
|
|
|
// 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.
|
|
|
|
//============================================================================
|
|
|
|
#ifndef vtk_m_count_ArrayHandleRandomStandardNormal_h
|
|
|
|
#define vtk_m_count_ArrayHandleRandomStandardNormal_h
|
|
|
|
|
|
|
|
#include <vtkm/Math.h>
|
|
|
|
#include <vtkm/cont/ArrayHandleRandomUniformReal.h>
|
|
|
|
#include <vtkm/cont/ArrayHandleTransform.h>
|
|
|
|
#include <vtkm/cont/ArrayHandleZip.h>
|
|
|
|
|
|
|
|
namespace vtkm
|
|
|
|
{
|
|
|
|
namespace cont
|
|
|
|
{
|
|
|
|
namespace detail
|
|
|
|
{
|
|
|
|
struct BoxMuller
|
|
|
|
{
|
|
|
|
VTKM_EXEC_CONT vtkm::Float32 operator()(const vtkm::Pair<vtkm::Float32, vtkm::Float32>& uv) const
|
|
|
|
{
|
|
|
|
// We take two U(0, 1) and return one N(0, 1)
|
|
|
|
return vtkm::Sqrt(-2.0f * vtkm::Log(uv.first)) * vtkm::Cos(2.0f * vtkm::TwoPif() * uv.second);
|
|
|
|
}
|
|
|
|
|
|
|
|
VTKM_EXEC_CONT vtkm::Float64 operator()(const vtkm::Pair<vtkm::Float64, vtkm::Float64>& uv) const
|
|
|
|
{
|
|
|
|
// We take two U(0, 1) and return one N(0, 1)
|
|
|
|
return vtkm::Sqrt(-2.0 * vtkm::Log(uv.first)) * vtkm::Cos(2 * vtkm::TwoPi() * uv.second);
|
|
|
|
}
|
|
|
|
};
|
|
|
|
} //detail
|
|
|
|
|
|
|
|
template <typename Real = vtkm::Float64>
|
|
|
|
class VTKM_ALWAYS_EXPORT ArrayHandleRandomStandardNormal
|
|
|
|
: public vtkm::cont::ArrayHandleTransform<
|
|
|
|
vtkm::cont::ArrayHandleZip<vtkm::cont::ArrayHandleRandomUniformReal<Real>,
|
|
|
|
vtkm::cont::ArrayHandleRandomUniformReal<Real>>,
|
|
|
|
detail::BoxMuller>
|
|
|
|
{
|
|
|
|
public:
|
|
|
|
using SeedType = vtkm::Vec<vtkm::UInt32, 1>;
|
|
|
|
using UniformReal = vtkm::cont::ArrayHandleRandomUniformReal<Real>;
|
|
|
|
|
|
|
|
VTKM_ARRAY_HANDLE_SUBCLASS(
|
|
|
|
ArrayHandleRandomStandardNormal,
|
|
|
|
(ArrayHandleRandomStandardNormal<Real>),
|
|
|
|
(vtkm::cont::ArrayHandleTransform<
|
|
|
|
vtkm::cont::ArrayHandleZip<vtkm::cont::ArrayHandleRandomUniformReal<Real>,
|
|
|
|
vtkm::cont::ArrayHandleRandomUniformReal<Real>>,
|
|
|
|
detail::BoxMuller>));
|
|
|
|
|
|
|
|
explicit ArrayHandleRandomStandardNormal(vtkm::Id length,
|
|
|
|
SeedType seed = { std::random_device{}() })
|
2020-06-25 15:51:06 +00:00
|
|
|
: Superclass(vtkm::cont::make_ArrayHandleZip(UniformReal{ length, seed },
|
|
|
|
UniformReal{ length, { ~seed[0] } }),
|
2020-06-17 21:48:09 +00:00
|
|
|
detail::BoxMuller{})
|
|
|
|
{
|
|
|
|
}
|
|
|
|
};
|
|
|
|
}
|
|
|
|
}
|
|
|
|
#endif // vtk_m_count_ArrayHandleRandomStandardNormal_h
|