Commit 271579a6 authored by Dave Moxey's avatar Dave Moxey
Browse files

Merge branch 'fix/whitenoise' into 'master'

Fix FillWhiteNoise

This MR fixes a bug in `Vmath::FillWhiteNoise` which causes it to alternate between only 2 values.

Instead of trying to pass a unique seed on each call of FillWhiteNoise, I decided to partly recover an old implementation where `seed` was an static variable ( therefore making consecutive calls part of a single larger random sequence). From what I understood, this had been changed to allow specifying different seeds for each MPI process. With this in mind, we can still pass an optional seed to replace the stored value.

An example of how to use this new logic is to use (the negative of) the rank as seed in the first call to FillWhiteNoise to decouple the parallel processes, and not passing any seed in the subsequent calls . Alternatively, we can still pass unique seeds on each call.

See merge request !718
parents 1c300139 e7cd5cd4
......@@ -34,6 +34,8 @@ v4.4.0
- Enabled MUMPS support in PETSc if a Fortran compiler was found and added 3D
support to the Helmholtz smoother used e.g. in FieldConverts C0Projection
module (!714)
- Fix bug in `Vmath::FillWhiteNoise` which caused `ForcingNoise` to have
a repeated pattern (!718)
- Fix bug in the calculation of the RHS magnitude in CG solver (!721)
- Fix bug in CMake Homebrew and MacPorts detection for OS X (!729)
- Fix bug in FieldUtils when using half mode expansions (!734)
......
......@@ -135,30 +135,44 @@ namespace Vmath
template LIB_UTILITIES_EXPORT Nektar::NekDouble ran2 (long* idum);
/// \brief Fills a vector with white noise.
template<class T> void FillWhiteNoise( int n, const T eps, T *x, const int incx, int outseed)
template<class T> void FillWhiteNoise( int n, const T eps, T *x,
const int incx, int outseed)
{
// Protect the static vars here and in ran2
boost::mutex::scoped_lock l(mutex);
// Define static variables for generating random numbers
static int iset = 0;
static T gset;
static long seed = 0;
// Bypass seed if outseed was specified
if( outseed != 9999)
{
seed = long(outseed);
}
while( n-- )
{
static int iset = 0;
static T gset;
long seed = long(outseed);
T fac, rsq, v1, v2;
if (iset == 0) {
do {
v1 = 2.0 * ran2<T> (&seed) - 1.0;
v2 = 2.0 * ran2<T> (&seed) - 1.0;
rsq = v1*v1 + v2*v2;
} while (rsq >= 1.0 || rsq == 0.0);
fac = sqrt(-2.0 * log (rsq) / rsq);
gset = v1 * fac;
iset = 1;
*x = eps * v2 * fac;
} else {
iset = 0;
*x = eps * gset;
if (iset == 0)
{
do
{
v1 = 2.0 * ran2<T> (&seed) - 1.0;
v2 = 2.0 * ran2<T> (&seed) - 1.0;
rsq = v1*v1 + v2*v2;
} while (rsq >= 1.0 || rsq == 0.0);
fac = sqrt(-2.0 * log (rsq) / rsq);
gset = v1 * fac;
iset = 1;
*x = eps * v2 * fac;
}
else
{
iset = 0;
*x = eps * gset;
}
x += incx;
}
......
......@@ -54,7 +54,8 @@ namespace Vmath
/// \brief Fills a vector with white noise.
template<class T> LIB_UTILITIES_EXPORT void FillWhiteNoise( int n, const T eps, T *x, const int incx, int seed = 0);
template<class T> LIB_UTILITIES_EXPORT void FillWhiteNoise(
int n, const T eps, T *x, const int incx, int seed = 9999);
/// \brief Multiply vector z = x*y
template<class T> LIB_UTILITIES_EXPORT void Vmul( int n, const T *x, const int incx, const T *y,
......
......@@ -53,7 +53,9 @@
Fill(n,alpha,&x[0],incx);
}
template<class T> void FillWhiteNoise( int n, const T eps, Array<OneD, T> &x, const int incx, int outseed)
template<class T> void FillWhiteNoise(
int n, const T eps, Array<OneD, T> &x,
const int incx, int outseed = 9999)
{
ASSERTL1(n*incx <= x.num_elements()+x.GetOffset(),"Out of bounds");
......
......@@ -97,11 +97,14 @@ namespace SolverUtils
m_Forcing = Array<OneD, Array<OneD, NekDouble> > (m_NumVariable);
// Fill forcing: use -i as seed to avoid repeated results
for (int i = 0; i < m_NumVariable; ++i)
// Fill forcing: use rank in seed to avoid repeated results
int seed = - m_session->GetComm()->GetRank();
m_Forcing[0] = Array<OneD, NekDouble> (nq, 0.0);
Vmath::FillWhiteNoise(nq,m_noise,m_Forcing[0],1,seed);
for (int i = 1; i < m_NumVariable; ++i)
{
m_Forcing[i] = Array<OneD, NekDouble> (nq, 0.0);
Vmath::FillWhiteNoise(nq,m_noise,&m_Forcing[i][0],1,-i);
Vmath::FillWhiteNoise(nq,m_noise,m_Forcing[i],1);
}
m_index = 0;
......@@ -125,8 +128,7 @@ namespace SolverUtils
for (int i = 0; i < m_NumVariable; ++i)
{
Vmath::FillWhiteNoise(outarray[i].num_elements(),
m_noise,&m_Forcing[i][0],1,
-m_index*m_NumVariable-i);
m_noise,m_Forcing[i],1);
}
}
......
......@@ -709,10 +709,12 @@ namespace Nektar
if (Noise > 0.0)
{
int seed = - m_comm->GetRank()*m_nConvectiveFields;
for (int i = 0; i < m_nConvectiveFields; i++)
{
Vmath::FillWhiteNoise(phystot, Noise, noise, 1,
m_comm->GetColumnComm()->GetRank()+1);
seed);
--seed;
Vmath::Vadd(phystot, m_fields[i]->GetPhys(), 1,
noise, 1, m_fields[i]->UpdatePhys(), 1);
m_fields[i]->FwdTrans_IterPerExp(m_fields[i]->GetPhys(),
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000483502</value>
<value variable="rho" tolerance="1e-12">0.000483148</value>
<value variable="rhou" tolerance="1e-12">48.1265</value>
<value variable="rhov" tolerance="1e-12">0.17737</value>
<value variable="rhow" tolerance="1e-12">7.15155e-06</value>
<value variable="rhov" tolerance="1e-12">0.177369</value>
<value variable="rhow" tolerance="1e-12">7.074e-06</value>
<value variable="E" tolerance="1e-12">17518.9</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00208788</value>
<value variable="rho" tolerance="1e-12">0.00207205</value>
<value variable="rhou" tolerance="1e-12">83.3318</value>
<value variable="rhov" tolerance="1e-12">0.752415</value>
<value variable="rhow" tolerance="1e-12">4.89968e-05</value>
<value variable="rhov" tolerance="1e-12">0.7524</value>
<value variable="rhow" tolerance="1e-12">3.34657e-05</value>
<value variable="E" tolerance="1e-12">18726.1</value>
</metric>
</metrics>
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000483502</value>
<value variable="rho" tolerance="1e-12">0.000483148</value>
<value variable="rhou" tolerance="1e-12">48.1265</value>
<value variable="rhov" tolerance="1e-12">0.17737</value>
<value variable="rhow" tolerance="1e-12"> 7.15155e-06</value>
<value variable="rhov" tolerance="1e-12">0.177369</value>
<value variable="rhow" tolerance="1e-12"> 7.074e-06</value>
<value variable="E" tolerance="1e-12">17518.9</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00208788</value>
<value variable="rho" tolerance="1e-12">0.00207205</value>
<value variable="rhou" tolerance="1e-12">83.3318</value>
<value variable="rhov" tolerance="1e-12">0.752415</value>
<value variable="rhow" tolerance="1e-12">4.89968e-05</value>
<value variable="rhov" tolerance="1e-12">0.7524</value>
<value variable="rhow" tolerance="1e-12">3.34657e-05</value>
<value variable="E" tolerance="1e-12">18726.1</value>
</metric>
</metrics>
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000397706</value>
<value variable="rho" tolerance="1e-12">0.000397131</value>
<value variable="rhou" tolerance="1e-12">48.1295</value>
<value variable="rhov" tolerance="1e-8">0.145836</value>
<value variable="rhow" tolerance="1e-8">8.79652e-06</value>
<value variable="rhov" tolerance="1e-8">0.145835</value>
<value variable="rhow" tolerance="1e-8">9.3349e-06</value>
<value variable="E" tolerance="1e-12">17519.9</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00139611</value>
<value variable="rhou" tolerance="1e-12">83.3516</value>
<value variable="rhov" tolerance="1e-8">0.505196</value>
<value variable="rhow" tolerance="1e-8">3.177e-05</value>
<value variable="rho" tolerance="1e-12">0.00139966</value>
<value variable="rhou" tolerance="1e-12">83.3517</value>
<value variable="rhov" tolerance="1e-8">0.50519</value>
<value variable="rhow" tolerance="1e-8">3.11563e-05</value>
<value variable="E" tolerance="1e-12">18953</value>
</metric>
</metrics>
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000397706</value>
<value variable="rho" tolerance="1e-12">0.000397131</value>
<value variable="rhou" tolerance="1e-12">48.1295</value>
<value variable="rhov" tolerance="1e-8">0.145836</value>
<value variable="rhow" tolerance="1e-8">8.79652e-06</value>
<value variable="rhov" tolerance="1e-8">0.145835</value>
<value variable="rhow" tolerance="1e-8">9.3349e-06</value>
<value variable="E" tolerance="1e-12">17519.9</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00139611</value>
<value variable="rhou" tolerance="1e-12">83.3516</value>
<value variable="rhov" tolerance="1e-8">0.505196</value>
<value variable="rhow" tolerance="1e-8">3.177e-05</value>
<value variable="rho" tolerance="1e-12">0.00139966</value>
<value variable="rhou" tolerance="1e-12">83.3517</value>
<value variable="rhov" tolerance="1e-8">0.50519</value>
<value variable="rhow" tolerance="1e-8">3.11563e-05</value>
<value variable="E" tolerance="1e-12">18953</value>
</metric>
</metrics>
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000391612</value>
<value variable="rho" tolerance="1e-12">0.000391074</value>
<value variable="rhou" tolerance="1e-12">48.1285</value>
<value variable="rhov" tolerance="1e-12">0.143446</value>
<value variable="rhow" tolerance="1e-12"> 6.98852e-06</value>
<value variable="rhov" tolerance="1e-12">0.143445</value>
<value variable="rhow" tolerance="1e-12"> 6.90409e-06</value>
<value variable="E" tolerance="1e-12">17519.5</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00164134</value>
<value variable="rho" tolerance="1e-12">0.00162551</value>
<value variable="rhou" tolerance="1e-12">83.3449</value>
<value variable="rhov" tolerance="1e-12">0.5882</value>
<value variable="rhow" tolerance="1e-12">4.62189e-05</value>
<value variable="rhov" tolerance="1e-12">0.588196</value>
<value variable="rhow" tolerance="1e-12">3.33264e-05</value>
<value variable="E" tolerance="1e-12">18878.3</value>
</metric>
</metrics>
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000391612</value>
<value variable="rho" tolerance="1e-12">0.000391074</value>
<value variable="rhou" tolerance="1e-12">48.1285</value>
<value variable="rhov" tolerance="1e-12">0.143446</value>
<value variable="rhow" tolerance="1e-12"> 6.98852e-06</value>
<value variable="rhov" tolerance="1e-12">0.143445</value>
<value variable="rhow" tolerance="1e-12"> 6.90409e-06</value>
<value variable="E" tolerance="1e-12">17519.5</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00164134</value>
<value variable="rho" tolerance="1e-12">0.00162551</value>
<value variable="rhou" tolerance="1e-12">83.3449</value>
<value variable="rhov" tolerance="1e-12">0.5882</value>
<value variable="rhow" tolerance="1e-12">4.62189e-05</value>
<value variable="rhov" tolerance="1e-12">0.588196</value>
<value variable="rhow" tolerance="1e-12">3.33264e-05</value>
<value variable="E" tolerance="1e-12">18878.3</value>
</metric>
</metrics>
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000397706</value>
<value variable="rho" tolerance="1e-12">0.000397131</value>
<value variable="rhou" tolerance="1e-12">48.1295</value>
<value variable="rhov" tolerance="1e-12">0.145836</value>
<value variable="rhow" tolerance="1e-8">8.79652e-06</value>
<value variable="rhov" tolerance="1e-12">0.145835</value>
<value variable="rhow" tolerance="1e-8">9.33488e-06</value>
<value variable="E" tolerance="1e-12">17519.9</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00139611</value>
<value variable="rhou" tolerance="1e-12">83.3516</value>
<value variable="rhov" tolerance="1e-12">0.505196</value>
<value variable="rhow" tolerance="1e-8">3.177e-05</value>
<value variable="rho" tolerance="1e-12">0.00139966</value>
<value variable="rhou" tolerance="1e-12">83.3517</value>
<value variable="rhov" tolerance="1e-12">0.50519</value>
<value variable="rhow" tolerance="1e-8">3.11565e-05</value>
<value variable="E" tolerance="1e-12">18953</value>
</metric>
</metrics>
......
......@@ -8,17 +8,17 @@
</files>
<metrics>
<metric type="L2" id="1">
<value variable="rho" tolerance="1e-12">0.000397706</value>
<value variable="rho" tolerance="1e-12">0.000397131</value>
<value variable="rhou" tolerance="1e-12">48.1295</value>
<value variable="rhov" tolerance="1e-12">0.145836</value>
<value variable="rhow" tolerance="1e-12">8.79677e-06</value>
<value variable="rhov" tolerance="1e-12">0.145835</value>
<value variable="rhow" tolerance="1e-12">9.33488e-06</value>
<value variable="E" tolerance="1e-12">17519.9</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.00139611</value>
<value variable="rhou" tolerance="1e-12">83.3516</value>
<value variable="rhov" tolerance="1e-12">0.505196</value>
<value variable="rhow" tolerance="1e-12">3.17716e-05</value>
<value variable="rho" tolerance="1e-12">0.00139966</value>
<value variable="rhou" tolerance="1e-12">83.3517</value>
<value variable="rhov" tolerance="1e-12">0.50519</value>
<value variable="rhow" tolerance="1e-12">3.11565e-05</value>
<value variable="E" tolerance="1e-12">18953</value>
</metric>
</metrics>
......
......@@ -15,11 +15,11 @@
<value variable="E" tolerance="1e-12">4.43372e+06</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.297684</value>
<value variable="rho" tolerance="1e-12">0.297688</value>
<value variable="rhou" tolerance="1e-12">86.8671</value>
<value variable="rhov" tolerance="1e-12"> 15.3684</value>
<value variable="rhow" tolerance="1e-12">1.00005</value>
<value variable="E" tolerance="1e-12">277512</value>
<value variable="rhow" tolerance="1e-12">1.00004</value>
<value variable="E" tolerance="1e-12">277513</value>
</metric>
</metrics>
</test>
......
......@@ -15,11 +15,11 @@
<value variable="E" tolerance="1e-12">4.43372e+06</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.297684</value>
<value variable="rho" tolerance="1e-12">0.297688</value>
<value variable="rhou" tolerance="1e-12">86.8671</value>
<value variable="rhov" tolerance="1e-12"> 15.3684</value>
<value variable="rhow" tolerance="1e-12">1.00005</value>
<value variable="E" tolerance="1e-12">277512</value>
<value variable="rhow" tolerance="1e-12">1.00004</value>
<value variable="E" tolerance="1e-12">277513</value>
</metric>
</metrics>
</test>
......
......@@ -15,7 +15,7 @@
<value variable="E" tolerance="1e-12">4.43372e+06</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12"> 0.275032</value>
<value variable="rho" tolerance="1e-12"> 0.275036</value>
<value variable="rhou" tolerance="1e-12">82.5652</value>
<value variable="rhov" tolerance="1e-8">10.8489</value>
<value variable="rhow" tolerance="1e-12">1.00004</value>
......
......@@ -15,7 +15,7 @@
<value variable="E" tolerance="1e-12">4.43372e+06</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12"> 0.275032</value>
<value variable="rho" tolerance="1e-12"> 0.275036</value>
<value variable="rhou" tolerance="1e-12">82.5652</value>
<value variable="rhov" tolerance="1e-8">10.8489</value>
<value variable="rhow" tolerance="1e-12">1.00004</value>
......
......@@ -15,10 +15,10 @@
<value variable="E" tolerance="1e-12">4.43372e+06</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.281922</value>
<value variable="rho" tolerance="1e-12">0.281925</value>
<value variable="rhou" tolerance="1e-12">85.4194</value>
<value variable="rhov" tolerance="1e-8">12.2322</value>
<value variable="rhow" tolerance="1e-12">1.00005</value>
<value variable="rhow" tolerance="1e-12">1.00004</value>
<value variable="E" tolerance="1e-12">272860</value>
</metric>
</metrics>
......
......@@ -15,10 +15,10 @@
<value variable="E" tolerance="1e-12">4.43372e+06</value>
</metric>
<metric type="Linf" id="2">
<value variable="rho" tolerance="1e-12">0.281922</value>
<value variable="rho" tolerance="1e-12">0.281925</value>
<value variable="rhou" tolerance="1e-12">85.4194</value>
<value variable="rhov" tolerance="1e-8">12.2322</value>
<value variable="rhow" tolerance="1e-12">1.00005</value>
<value variable="rhow" tolerance="1e-12">1.00004</value>
<value variable="E" tolerance="1e-12">272860</value>
</metric>
</metrics>
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment