1 文本格式
using System;
using System.Collections.Generic;
namespace Legalsoft.Truffer
{
/// <summary>
/// Sobol quasi-random sequence
/// </summary>
public class Sobol
{
public Sobol() { }
public static void sobseq(int n, double[] x)
{
const int MAXBIT = 30;
const int MAXDIM = 6;
uint xin = 0;
int[] mdeg = new int[] { 1, 2, 3, 3, 4, 4 };
uint[] ix = new uint[MAXDIM];
List<uint[]> iu = new List<uint[]>(MAXBIT);
uint[] ip = new uint[] {
0, 1, 1, 2, 1, 4
};
uint[] iv = new uint[] {
1, 1, 1, 1, 1, 1,
3, 1, 3, 3, 1, 1,
5, 7, 7, 3, 3, 5,
15, 11, 5, 15, 13, 9
};
double fac = 0.0;
if (n < 0)
{
for (int k = 0; k < MAXDIM; k++)
{
ix[k] = 0;
}
xin = 0;
if (iv[0] != 1)
{
return;
}
fac = 1.0 / (1 << MAXBIT);
for (int j = 0, k = 0; j < MAXBIT; j++, k += MAXDIM)
{
iu[j] = new uint[iv[k]];
}
for (int k = 0; k < MAXDIM; k++)
{
for (int j = 0; j < mdeg[k]; j++)
{
iu[j][k] <<= (MAXBIT - 1 - j);
}
for (int j = mdeg[k]; j < MAXBIT; j++)
{
uint ipp = ip[k];
uint i = iu[j - mdeg[k]][k];
i ^= (i >> mdeg[k]);
for (int l = mdeg[k] - 1; l >= 1; l--)
{
if ((ipp & 1) != 0)
{
i ^= iu[j - l][k];
}
ipp >>= 1;
}
iu[j][k] = i;
}
}
}
else
{
uint im = xin++;
int j = 0;
for (; j < MAXBIT; j++)
{
if ((im & 1) == 0)
{
break;
}
im >>= 1;
}
if (j >= MAXBIT)
{
throw new Exception("MAXBIT too small in sobseq");
}
im = (uint)(j * MAXDIM);
for (int k = 0; k < Math.Min(n, MAXDIM); k++)
{
ix[k] ^= iv[im + k];
x[k] = ix[k] * fac;
}
}
}
}
}
2 代码格式
using System;
using System.Collections.Generic;
namespace Legalsoft.Truffer
{
/// <summary>
/// Sobol quasi-random sequence
/// </summary>
public class Sobol
{
public Sobol() { }
public static void sobseq(int n, double[] x)
{
const int MAXBIT = 30;
const int MAXDIM = 6;
uint xin = 0;
int[] mdeg = new int[] { 1, 2, 3, 3, 4, 4 };
uint[] ix = new uint[MAXDIM];
List<uint[]> iu = new List<uint[]>(MAXBIT);
uint[] ip = new uint[] {
0, 1, 1, 2, 1, 4
};
uint[] iv = new uint[] {
1, 1, 1, 1, 1, 1,
3, 1, 3, 3, 1, 1,
5, 7, 7, 3, 3, 5,
15, 11, 5, 15, 13, 9
};
double fac = 0.0;
if (n < 0)
{
for (int k = 0; k < MAXDIM; k++)
{
ix[k] = 0;
}
xin = 0;
if (iv[0] != 1)
{
return;
}
fac = 1.0 / (1 << MAXBIT);
for (int j = 0, k = 0; j < MAXBIT; j++, k += MAXDIM)
{
iu[j] = new uint[iv[k]];
}
for (int k = 0; k < MAXDIM; k++)
{
for (int j = 0; j < mdeg[k]; j++)
{
iu[j][k] <<= (MAXBIT - 1 - j);
}
for (int j = mdeg[k]; j < MAXBIT; j++)
{
uint ipp = ip[k];
uint i = iu[j - mdeg[k]][k];
i ^= (i >> mdeg[k]);
for (int l = mdeg[k] - 1; l >= 1; l--)
{
if ((ipp & 1) != 0)
{
i ^= iu[j - l][k];
}
ipp >>= 1;
}
iu[j][k] = i;
}
}
}
else
{
uint im = xin++;
int j = 0;
for (; j < MAXBIT; j++)
{
if ((im & 1) == 0)
{
break;
}
im >>= 1;
}
if (j >= MAXBIT)
{
throw new Exception("MAXBIT too small in sobseq");
}
im = (uint)(j * MAXDIM);
for (int k = 0; k < Math.Min(n, MAXDIM); k++)
{
ix[k] ^= iv[im + k];
x[k] = ix[k] * fac;
}
}
}
}
}