diff --git a/example.conf b/example.conf index e5a7db4..e383bba 100644 --- a/example.conf +++ b/example.conf @@ -107,7 +107,8 @@ seed = 12345 ##> The MUSIC1 multi-scale random number generator is provided for convenience ## warning: MUSIC1 generator is not MPI parallel (yet) (memory is needed for full field on each task) -# generator = MUSIC1 +# generator = MUSIC1 +# music2_rng = false # seed[7] = 12345 # seed[8] = 23456 # seed[9] = 34567 diff --git a/src/plugins/random_music.cc b/src/plugins/random_music.cc index 675c318..4fdc8a4 100644 --- a/src/plugins/random_music.cc +++ b/src/plugins/random_music.cc @@ -216,6 +216,11 @@ void RNG_music::parse_random_parameters(void) void RNG_music::compute_random_numbers(void) { bool rndsign = pcf_->get_value_safe("random", "grafic_sign", false); + bool music2_rng = pcf_->get_value_safe("random", "music2_rng", false); + if (music2_rng) + { + music::ilog.Print("MUSIC1 RNG plugin: using MUSIC2 compatible seeds"); + } //--- FILL ALL WHITE NOISE ARRAYS FOR WHICH WE NEED THE FULL FIELD ---// @@ -235,9 +240,9 @@ void RNG_music::compute_random_numbers(void) //#warning add possibility to read noise from file also here! if (rngfnames_[i].size() > 0) - music::ilog.Print("Warning: Cannot use filenames for higher levels currently! Ignoring!"); + music::wlog.Print("Warning: Cannot use filenames for higher levels currently! Ignoring!"); - randc_[i] = new rng(*randc_[i - 1], ran_cube_size_, rngseeds_[i], true); + randc_[i] = new rng(*randc_[i - 1], ran_cube_size_, rngseeds_[i], music2_rng, true); delete randc_[i - 1]; randc_[i - 1] = NULL; } diff --git a/src/plugins/random_music_wnoise_generator.cc b/src/plugins/random_music_wnoise_generator.cc index 80d85d9..a0a169e 100644 --- a/src/plugins/random_music_wnoise_generator.cc +++ b/src/plugins/random_music_wnoise_generator.cc @@ -22,6 +22,21 @@ #include #include "random_music_wnoise_generator.hh" +inline double Meyer_scaling_function( double k, double kmax ) +{ + constexpr double twopithirds{2.0*M_PI/3.0}; + constexpr double fourpithirds{4.0*M_PI/3.0}; + auto nu = []( double x ){ return x<0.0?0.0:(x<1.0?x:1.0); }; + + k = std::abs(k)/kmax * 2 * M_PI; + + if( k < twopithirds ) return 1.0; + else if( k< fourpithirds ){ + return std::cos( 0.5*M_PI * nu(3*k/(2*M_PI)-1.0) ); + } + return 0.0; +} + template music_wnoise_generator::music_wnoise_generator(unsigned res, unsigned cubesize, long baseseed, int *x0, int *lx) : res_(res), cubesize_(cubesize), ncubes_(1), baseseed_(baseseed) @@ -43,7 +58,7 @@ music_wnoise_generator::music_wnoise_generator(unsigned res, unsigned cubesiz initialize(); mean = fill_all(); - + if (zeromean) { mean = 0.0; @@ -389,7 +404,7 @@ music_wnoise_generator::music_wnoise_generator(/*const*/ music_wnoise_generat template music_wnoise_generator::music_wnoise_generator(music_wnoise_generator &rc, unsigned cubesize, long baseseed, - bool kspace, bool isolated, int *x0_, int *lx_, bool zeromean) + bool music2_rng, bool kspace, bool isolated, int *x0_, int *lx_, bool zeromean) : res_(2 * rc.res_), cubesize_(cubesize), ncubes_(1), baseseed_(baseseed) { initialize(); @@ -496,15 +511,16 @@ music_wnoise_generator::music_wnoise_generator(music_wnoise_generator &rc, val *= val_phas * sqrt8; - if (i != (int)nxc / 2 && j != (int)nyc / 2 && k != (int)nzc / 2) - { - cfine[qf][0] = val.real(); - cfine[qf][1] = val.imag(); - } - else - { - //RE(cfine[qf]) = val.real(); - //IM(cfine[qf]) = 0.0; + if(i != (int)nxc / 2 && j != (int)nyc / 2 && k != (int)nzc / 2){ + double blend_coarse_x = Meyer_scaling_function(kx, nxc / 2); + double blend_coarse_y = Meyer_scaling_function(ky, nyc / 2); + double blend_coarse_z = Meyer_scaling_function(kz, nzc / 2); + + double blend_coarse = blend_coarse_x*blend_coarse_y*blend_coarse_z; + double blend_fine = std::sqrt(1.0-blend_coarse*blend_coarse); + + cfine[qf][0] = blend_fine * cfine[qf][0] + blend_coarse * val.real(); + cfine[qf][1] = blend_fine * cfine[qf][1] + blend_coarse * val.imag(); } } diff --git a/src/plugins/random_music_wnoise_generator.hh b/src/plugins/random_music_wnoise_generator.hh index 0ab3a07..a39ef7f 100644 --- a/src/plugins/random_music_wnoise_generator.hh +++ b/src/plugins/random_music_wnoise_generator.hh @@ -162,7 +162,7 @@ public: music_wnoise_generator(unsigned res, unsigned cubesize, long baseseed, int *x0, int *lx); //! constructor for constrained fine field - music_wnoise_generator(music_wnoise_generator &rc, unsigned cubesize, long baseseed, + music_wnoise_generator(music_wnoise_generator &rc, unsigned cubesize, long baseseed, bool music2_rng, bool kspace = false, bool isolated = false, int *x0_ = NULL, int *lx_ = NULL, bool zeromean = true); //! constructor