wip

WeiqunZhang · Aug 12, 2023 · 7a7e4ae · 7a7e4ae
1 parent 9bceb99
commit 7a7e4ae
Show file tree

Hide file tree

Showing 6 changed files with 115 additions and 1 deletion.
diff --git a/Src/Base/AMReX_GpuError.H b/Src/Base/AMReX_GpuError.H
@@ -79,6 +79,11 @@ namespace Gpu {
  + ": " + cudaGetErrorString(amrex_i_err)); \
  amrex::Abort(errStr); \
  }}
+
+#define AMREX_CURAND_SAFE_CALL(x) do { if((x)!=CURAND_STATUS_SUCCESS) { \
+ std::string errStr(std::string("CURAND error in file ") + __FILE__ \
+ + " line " + std::to_string(__LINE__)); \
+ amrex::Abort(errStr); }} while(0)
 #endif
 
 #ifdef AMREX_USE_HIP

diff --git a/Src/Base/AMReX_MultiFabUtil.H b/Src/Base/AMReX_MultiFabUtil.H
@@ -361,6 +361,31 @@ namespace amrex
  void FourthOrderInterpFromFineToCoarse (MultiFab& cmf, int scomp, int ncomp,
  MultiFab const& fmf,
  IntVect const& ratio);
+
+ /**
+ * \brief Fill MultiFab with random numbers from uniform distribution
+ *
+ * The uniform distribution range is [0.0, 1.0). All cells including
+ * ghost cells are filled.
+ *
+ * \param mf MultiFab
+ * \param scomp starting component
+ * \param ncomp number of component
+ */
+ void FillRandom (MultiFab& mf, int scomp, int ncomp);
+
+ /**
+ * \brief Fill MultiFab with random numbers from nornmal distribution
+ *
+ * All cells including ghost cells are filled.
+ *
+ * \param mf MultiFab
+ * \param scomp starting component
+ * \param ncomp number of component
+ * \param mean mean of normal distribution
+ * \param stddev standard deviation of normal distribution
+ */
+ void FillRandomNormal (MultiFab& mf, int scomp, int ncomp, Real mean, Real stddev);
 }
 
 namespace amrex {

diff --git a/Src/Base/AMReX_MultiFabUtil.cpp b/Src/Base/AMReX_MultiFabUtil.cpp
@@ -1,5 +1,6 @@
 
 #include <AMReX_MultiFabUtil.H>
+#include <AMReX_Random.H>
 #include <sstream>
 #include <iostream>
 
@@ -1184,4 +1185,30 @@ namespace amrex
 
  cmf.ParallelCopy(tmp, 0, scomp, ncomp);
  }
+
+ void FillRandom (MultiFab& mf, int scomp, int ncomp)
+ {
+#ifdef AMREX_USE_OMP
+#pragma omp parallel if (Gpu::notInLaunchRegion())
+#endif
+ for (MFIter mfi(mf); mfi.isValid(); ++mfi)
+ {
+ auto* p = mf[mfi].dataPtr(scomp);
+ Long npts = mf[mfi].box().numPts() * ncomp;
+ FillRandom(p, npts);
+ }
+ }
+
+ void FillRandomNormal (MultiFab& mf, int scomp, int ncomp, Real mean, Real stddev)
+ {
+#ifdef AMREX_USE_OMP
+#pragma omp parallel if (Gpu::notInLaunchRegion())
+#endif
+ for (MFIter mfi(mf); mfi.isValid(); ++mfi)
+ {
+ auto* p = mf[mfi].dataPtr(scomp);
+ Long npts = mf[mfi].box().numPts() * ncomp;
+ FillRandomNormal(p, npts, mean, stddev);
+ }
+ }
 }
diff --git a/Src/Base/AMReX_Random.H b/Src/Base/AMReX_Random.H
@@ -144,6 +144,12 @@ namespace amrex
  */
  ULong Random_long (ULong n); // [0,n-1]
 
+ //! Fill random numbers from uniform distribution
+ void FillRandom (Real* p, Long N);
+
+ //! Fill random numbers from normal distribution
+ void FillRandomNormal (Real* p, Long N, Real mean, Real stddev);
+
  namespace detail {
  inline ULong DefaultGpuSeed () {
  return ParallelDescriptor::MyProc()*1234567ULL + 12345ULL;

diff --git a/Src/Base/AMReX_Random.cpp b/Src/Base/AMReX_Random.cpp
@@ -25,6 +25,16 @@ namespace amrex {
 }
 #endif
 
+#ifdef AMREX_USE_GPU
+namespace
+{
+ bool generator_initialized = false;
+#if defined(AMREX_USE_CUDA)
+ curandGenerator_t cuda_rand_gen;
+#endif
+}
+#endif
+
 #ifdef AMREX_USE_GPU
 namespace {
 void ResizeRandomSeed (amrex::ULong gpu_seed)
@@ -204,6 +214,47 @@ DeallocateRandomSeedDevArray ()
 #endif
 }
 
+void FillRandom (Real* p, Long N)
+{
+#ifdef AMREX_USE_GPU
+#else
+ std::uniform_real_distribution<Real> distribution(Real(0.0), Real(1.0));
+ auto& gen = generators[OpenMP::get_thread_num()];
+ for (Long i = 0; i < N; ++i) {
+ p[i] = distribution(gen);
+ }
+#endif
+}
+
+void FillRandomNormal (Real* p, Long N, Real mean, Real stddev)
+{
+#ifdef AMREX_USE_GPU
+
+#if defined(AMREX_USE_CUDA)
+ if (! generator_initialized) {
+ AMREX_CURAND_SAFE_CALL(curandCreateGenerator(&cuda_rand_gen,
+ CURAND_RNG_PSEUDO_DEFAULT));
+ AMREX_CURAND_SAFE_CALL(curandSetPseudoRandomGeneratorSeed(cuda_rand_gen,
+ 1234ULL));
+ }
+#ifdef BL_USE_FLOAT
+ AMREX_CURAND_SAFE_CALL(curandGenerateNormal(cuda_rand_gen, p, N, mean, stddev));
+#else
+ AMREX_CURAND_SAFE_CALL(curandGenerateNormalDouble(cuda_rand_gen, p, N, mean, stddev));
+#endif
+#endif
+
+ Gpu::streamSynchronize();
+
+#else
+ std::normal_distribution<Real> distribution(mean, stddev);
+ auto& gen = generators[OpenMP::get_thread_num()];
+ for (Long i = 0; i < N; ++i) {
+ p[i] = distribution(gen);
+ }
+#endif
+}
+
 } // namespace amrex
 
 

diff --git a/Tools/GNUMake/Make.defs b/Tools/GNUMake/Make.defs
@@ -1113,7 +1113,7 @@ else ifeq ($(USE_CUDA),TRUE)
  MAKE_CUDA_PATH := $(SYSTEM_CUDA_PATH)
  endif
 
- LIBRARIES += -lcuda
+ LIBRARIES += -lcuda -lcurand
 
  ifneq ($(MAKE_CUDA_PATH),)
  LIBRARY_LOCATIONS += $(MAKE_CUDA_PATH)/lib64