waves
Basic FE playground
wSpectralApproach.h
Go to the documentation of this file.
1 #ifndef KATOPTRON_SPECTRALAPPROACH_H
2 #define KATOPTRON_SPECTRALAPPROACH_H
3 
4 #include "EnsembleTraits.h"
5 #include <cmath>
6 #include <random>
7 
13 template <typename Scalar, typename Device>
14 Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> construct_Z(size_t muw, unsigned long seed_Z, size_t wait = 0)
15 {
16  Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> Z = Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device>("Z", muw, muw);
17 
18  std::default_random_engine generator;
19  generator.seed(seed_Z);
20  std::uniform_real_distribution<double> distribution(0.0, 1.0);
21  // Redraw the data related to the previous samples to start at the good place
22  for (size_t i = 0; i < wait; ++i)
23  distribution(generator);
24 
25  for (size_t ell = 0; ell < EnsembleTraits<Scalar>::size; ++ell)
26  for (size_t i = 0; i < muw; ++i)
27  for (size_t j = 0; j < muw; ++j)
28  EnsembleTraits<Scalar>::coeff(Z(i, j), ell) = sqrt(-log(distribution(generator)));
29 
30  return Z;
31 }
32 
38 template <typename Scalar, typename Device>
39 Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> construct_Phi(size_t muw, unsigned long seed_Phi, size_t wait = 0)
40 {
41  Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> Phi = Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device>("Phi", muw, muw);
42 
43  std::default_random_engine generator;
44  generator.seed(seed_Phi);
45  std::uniform_real_distribution<double> distribution(0.0, 1.0);
46  // Redraw the data related to the previous samples to start at the good place
47  for (size_t i = 0; i < wait; ++i)
48  distribution(generator);
49 
50  for (size_t ell = 0; ell < EnsembleTraits<Scalar>::size; ++ell)
51  for (size_t i = 0; i < muw; ++i)
52  for (size_t j = 0; j < muw; ++j)
53  EnsembleTraits<Scalar>::coeff(Phi(i, j), ell) = 2. * M_PI * distribution(generator);
54 
55  return Phi;
56 }
57 
62 template <typename Scalar, typename Device>
64 {
65 public:
66  typedef typename Teuchos::ScalarTraits<Scalar>::coordinateType MeshScalar;
67 
68  Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> Z; // Z=sqrt(-log(rand(muw,muw)));
69  Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> Phi; // Phi=2*pi*rand(muw,muw);
70  Kokkos::View<MeshScalar *, Kokkos::LayoutLeft, Device> w;
71  Kokkos::View<MeshScalar *, Kokkos::LayoutLeft, Device> S;
72 
74  size_t muw;
75 
76  SpectralApproach(Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> _Z,
77  Kokkos::View<Scalar **, Kokkos::LayoutLeft, Device> _Phi,
78  MeshScalar _a = 20, MeshScalar _r = 8, size_t _muw = 512,
79  bool Gaussian = true) : Z(_Z), Phi(_Phi), a(_a), r(_r), muw(_muw)
80  {
81  MeshScalar wL = M_PI / a * r;
82  Dw = 2 * wL / muw;
83 
84  w = Kokkos::View<MeshScalar *, Kokkos::LayoutLeft, Device>("w", muw);
85  S = Kokkos::View<MeshScalar *, Kokkos::LayoutLeft, Device>("w", muw);
86 
87  for (size_t i = 0; i < muw; ++i)
88  w(i) = -wL + (i + 0.5) * Dw;
89 
90  if (Gaussian)
91  for (size_t i = 0; i < muw; ++i)
92  S(i) = a / sqrt(2. * M_PI) * exp(-pow(a * w(i), 2) / 2.);
93  else
94  for (size_t i = 0; i < muw; ++i)
95  S(i) = a / (M_PI * (1 + pow(a, 2) * pow(w(i), 2)));
96  }
97 
98  KOKKOS_INLINE_FUNCTION
99  Scalar operator()(const MeshScalar x, const MeshScalar y) const
100  {
101  //Teuchos::TimeMonitor LocalTimer (*Random_field_eval);
102 
103  Scalar output = 0.;
104 
105  for (size_t i = 0; i < muw; ++i)
106  for (size_t j = 0; j < muw; ++j)
107  output += sqrt(S(i) * S(j)) * Z(i, j) * cos(Phi(i, j) + w(i) * x + w(j) * y);
108 
109  output *= sqrt(2) * Dw;
110 
111  return output;
112  }
113 };
114 
115 #endif //KATOPTRON_SPECTRALAPPROACH_H
construct_Z
Kokkos::View< Scalar **, Kokkos::LayoutLeft, Device > construct_Z(size_t muw, unsigned long seed_Z, size_t wait=0)
A function which computes a random matrix starting at a given seed.
Definition: wSpectralApproach.h:14
SpectralApproach::a
MeshScalar a
Definition: wSpectralApproach.h:73
SpectralApproach::operator()
KOKKOS_INLINE_FUNCTION Scalar operator()(const MeshScalar x, const MeshScalar y) const
Definition: wSpectralApproach.h:99
SpectralApproach
A class that constructs realization(s) of Gaussian or exponential scalar random field based on the sp...
Definition: wSpectralApproach.h:63
SpectralApproach::Dw
MeshScalar Dw
Definition: wSpectralApproach.h:73
SpectralApproach::MeshScalar
Teuchos::ScalarTraits< Scalar >::coordinateType MeshScalar
Definition: wSpectralApproach.h:66
SpectralApproach::Z
Kokkos::View< Scalar **, Kokkos::LayoutLeft, Device > Z
Definition: wSpectralApproach.h:68
EnsembleTraits
Definition: EnsembleTraits.h:8
SpectralApproach::SpectralApproach
SpectralApproach(Kokkos::View< Scalar **, Kokkos::LayoutLeft, Device > _Z, Kokkos::View< Scalar **, Kokkos::LayoutLeft, Device > _Phi, MeshScalar _a=20, MeshScalar _r=8, size_t _muw=512, bool Gaussian=true)
Definition: wSpectralApproach.h:76
SpectralApproach::muw
size_t muw
Definition: wSpectralApproach.h:74
SpectralApproach::w
Kokkos::View< MeshScalar *, Kokkos::LayoutLeft, Device > w
Definition: wSpectralApproach.h:70
EnsembleTraits.h
SpectralApproach::Phi
Kokkos::View< Scalar **, Kokkos::LayoutLeft, Device > Phi
Definition: wSpectralApproach.h:69
SpectralApproach::S
Kokkos::View< MeshScalar *, Kokkos::LayoutLeft, Device > S
Definition: wSpectralApproach.h:71
SpectralApproach::r
MeshScalar r
Definition: wSpectralApproach.h:73
construct_Phi
Kokkos::View< Scalar **, Kokkos::LayoutLeft, Device > construct_Phi(size_t muw, unsigned long seed_Phi, size_t wait=0)
A function which computes a random matrix starting at a given seed.
Definition: wSpectralApproach.h:39