如何使用multithreadingC＃实现Sierat Of Eratosthenes？
我正在尝试使用Mutithreading实现Sieve Of Eratosthenes。 这是我的实施：
using System; using System.Collections.Generic; using System.Threading; namespace Sieve_Of_Eratosthenes { class Controller { public static int upperLimit = 1000000000; public static bool[] primeArray = new bool[upperLimit]; static void Main(string[] args) { DateTime startTime = DateTime.Now; Initialize initial1 = new Initialize(0, 249999999); Initialize initial2 = new Initialize(250000000, 499999999); Initialize initial3 = new Initialize(500000000, 749999999); Initialize initial4 = new Initialize(750000000, 999999999); initial1.thread.Join(); initial2.thread.Join(); initial3.thread.Join(); initial4.thread.Join(); int sqrtLimit = (int)Math.Sqrt(upperLimit); Sieve sieve1 = new Sieve(249999999); Sieve sieve2 = new Sieve(499999999); Sieve sieve3 = new Sieve(749999999); Sieve sieve4 = new Sieve(999999999); for (int i = 3; i < sqrtLimit; i += 2) { if (primeArray[i] == true) { int squareI = i * i; if (squareI 249999999 & squareI 499999999 & squareI 749999999 & squareI <= 999999999) { sieve4.set(i); sieve4.thread.Join(); } } } int count = 0; primeArray[2] = true; for (int i = 2; i < upperLimit; i++) { if (primeArray[i]) { count++; } } Console.WriteLine("Total: " + count); DateTime endTime = DateTime.Now; TimeSpan elapsedTime = endTime - startTime; Console.WriteLine("Elapsed time: " + elapsedTime.Seconds); } public class Initialize { public Thread thread; private int lowerLimit; private int upperLimit; public Initialize(int lowerLimit, int upperLimit) { this.lowerLimit = lowerLimit; this.upperLimit = upperLimit; thread = new Thread(this.InitializeArray); thread.Priority = ThreadPriority.Highest; thread.Start(); } private void InitializeArray() { for (int i = this.lowerLimit; i <= this.upperLimit; i++) { if (i % 2 == 0) { Controller.primeArray[i] = false; } else { Controller.primeArray[i] = true; } } } } public class Sieve { public Thread thread; public int i; private int upperLimit; public Sieve(int upperLimit) { this.upperLimit = upperLimit; } public void set(int i) { this.i = i; thread = new Thread(this.primeGen); thread.Start(); } public void primeGen() { for (int j = this.i * this.i; j <= this.upperLimit; j += i) { Controller.primeArray[j] = false; } } } } } 这需要30秒才能产生输出，有什么方法可以加快速度吗？
编辑：这是TPL实现：
public LinkedList GetPrimeList(int limit) { LinkedList primeList = new LinkedList(); bool[] primeArray = new bool[limit]; Console.WriteLine("Initialization started..."); Parallel.For(0, limit, i => { if (i % 2 == 0) { primeArray[i] = false; } else { primeArray[i] = true; } } ); Console.WriteLine("Initialization finished..."); /*for (int i = 0; i { lock (this) { if (primeArray[i]) { for (int j = i * i; j < limit; j += i) { primeArray[j] = false; } } } } ); Console.WriteLine("Operation finished..."); /*for (int i = 3; i < sqrtLimit; i += 2) { if (primeArray[i]) { for (int j = i * i; j { lock (this) { if (primeArray[i]) { //primeList.AddLast(i); count++; } } } ); Console.WriteLine("Counting finished..."); Console.WriteLine(count); /*for (int i = 3; i < limit; i++) { if (primeArray[i]) { primeList.AddLast(i); } }*/ return primeList; }
谢谢。
编辑：
我对这个问题的回答是：是的，你绝对可以使用任务并行库（TPL）来快速找到10亿个素数。 问题中的给定代码很慢，因为它没有有效地使用内存或多处理，并且最终输出也没有效率。
因此， 除了多处理之外 ，您还可以做很多事情来加速Eratosthenese的筛选，如下所示：
- 你筛选所有数字，偶数和奇数，它们都使用更多的内存（10亿个字节，10亿的范围），并且由于不必要的处理而变慢。 只使用两个是唯一偶数素数的事实，这样使数组只表示奇数素数，这将占存储器需求的一半，并将复合数量剔除操作的数量减少两倍以上，因此操作可能需要20秒你的机器的质量达到十亿。
- 复合数量剔除如此庞大的内存数组的部分原因是它大大超过了CPU缓存大小，因此许多内存访问以某种随机的方式进入主内存，这意味着剔除给定的复合数字表示可能需要超过一百个CPU时钟周期，而如果它们都在L1高速缓存中，它只需要一个周期，而在L2高速缓存中只需要大约四个周期; 并非所有访问都采用最坏的情况，但这肯定会减慢处理速度。 使用位打包数组来表示主要候选者将减少8倍的内存使用，并使最坏情况的访问不那么常见。 虽然访问单个位会有计算开销，但您会发现存在净增益，因为减少平均内存访问时间所节省的时间将大于此成本。 实现它的简单方法是使用BitArray而不是bool数组。 使用shift和“and”操作编写自己的位访问将比使用BitArray类更有效。 您会发现使用BitArray进行轻微保存，另外两个因素是使用此更改进行单个线程性能（可能大约十或十二秒）进行自己的位操作。
- 您找到的素数计数的输出效率不高，因为它需要数组访问和每个候选素数的if条件。 一旦你将筛选缓冲区作为一个数组打包的字位数组，你可以通过计数查找表（LUT）更有效地完成这项工作，它可以消除if条件， 每个位打包字只需要两次数组访问。 这样做，与剔除复合数字的时间相比，计算时间成为工作中可忽略不计的部分，为了进一步节省，可以将素数计数减少到8秒，达到10亿。
- 进一步减少处理过的主要候选者的数量可以是应用车轮分解的结果，其中从处理中去除了素数2,3和5的因子并且通过调整比特填充的方法也可以增加有效范围。给定大小的位缓冲区大约是另一个大约两倍。 尽管以进一步的计算复杂性为代价，这可以将复合数量剔除操作的数量减少多达三倍以上的另一个重要因素。
- 为了进一步减少内存使用，使内存访问更加高效，并为每页段的多处理做好准备，可以将工作划分为不大于L1或L2高速缓存大小的页面。 这要求将所有素数的基本素数表保持到最大质数候选者的平方根，并重新计算在给定页面段中剔除时使用的每个基本素数的起始地址参数，但这仍然更有效而不是使用巨大的剔除arrays。 实现此页面分段的另一个好处是，之后不必提前指定上筛分限制，而是可以在处理更多上页时根据需要扩展基本素数。 通过所有这些优化，您可以在大约2.5秒内产生高达10亿的素数。
- 最后，可以使用TPL或Threads对页面段进行多处理的最后一步，使用大约L2缓存大小（每个核心）的缓冲区大小将在双核非Hyper Threaded上产生因子2的额外增益（HT）老款处理器作为英特尔E7500 Core2Duo的一个执行时间，找到的素数为10亿左右，约为1.25秒左右。
我已经实现了一个multithreading的Eratosthenes筛子作为另一个线程的答案，以显示Atkin筛选对于Eratosthenes筛子没有任何优势。 它使用Task和TaskFactory中的任务并行库（TPL），因此至少需要DotNet Framework 4.我使用上面讨论的所有优化作为同一问题的替代答案进一步调整了该代码。 我在这里重新发布了经过调整的代码，添加了注释和更易于阅读的格式，如下所示：
using System; using System.Collections; using System.Collections.Generic; using System.Linq; using System.Threading; using System.Threading.Tasks; class UltimatePrimesSoE : IEnumerable { #region const and static readonly field's, private struct's and classes //one can get single threaded performance by setting NUMPRCSPCS = 1 static readonly uint NUMPRCSPCS = (uint)Environment.ProcessorCount + 1; //the L1CACHEPOW can not be less than 14 and is usually the two raised to the power of the L1 or L2 cache const int L1CACHEPOW = 14, L1CACHESZ = (1 << L1CACHEPOW), MXPGSZ = L1CACHESZ / 2; //for buffer ushort[] const uint CHNKSZ = 17; //this times BWHLWRDS (below) times two should not be bigger than the L2 cache in bytes //the 2,3,57 factorial wheel increment pattern, (sum) 48 elements long, starting at prime 19 position static readonly byte[] WHLPTRN = { 2,3,1,3,2,1,2,3,3,1,3,2,1,3,2,3,4,2,1,2,1,2,4,3, 2,3,1,2,3,1,3,3,2,1,2,3,1,3,2,1,2,1,5,1,5,1,2,1 }; const uint FSTCP = 11; static readonly byte[] WHLPOS; static readonly byte[] WHLNDX; //look up wheel position from index and vice versa static readonly byte[] WHLRNDUP; //to look up wheel rounded up index positon values, allow for overflow in size static readonly uint WCRC = WHLPTRN.Aggregate(0u, (acc, n) => acc + n); //small wheel circumference for odd numbers static readonly uint WHTS = (uint)WHLPTRN.Length; static readonly uint WPC = WHTS >> 4; //number of wheel candidates static readonly byte[] BWHLPRMS = { 2,3,5,7,11,13,17 }; const uint FSTBP = 19; //big wheel primes, following prime //the big wheel circumference expressed in number of 16 bit words as in a minimum bit buffer size static readonly uint BWHLWRDS = BWHLPRMS.Aggregate(1u, (acc, p) => acc * p) / 2 / WCRC * WHTS / 16; //page size and range as developed from the above static readonly uint PGSZ = MXPGSZ / BWHLWRDS * BWHLWRDS; static readonly uint PGRNG = PGSZ * 16 / WHTS * WCRC; //buffer size (multiple chunks) as produced from the above static readonly uint BFSZ = CHNKSZ * PGSZ, BFRNG = CHNKSZ * PGRNG; //number of uints even number of caches in chunk static readonly ushort[] MCPY; //a Master Copy page used to hold the lower base primes preculled version of the page struct Wst { public ushort msk; public byte mlt; public byte xtr; public ushort nxt; } static readonly byte[] PRLUT; /*Wheel Index Look Up Table */ static readonly Wst[] WSLUT; //Wheel State Look Up Table static readonly byte[] CLUT; // a Counting Look Up Table for very fast counting of primes class Bpa { //very efficient auto-resizing thread-safe read-only indexer class to hold the base primes array byte[] sa = new byte[0]; uint lwi = 0, lpd = 0; object lck = new object(); public uint this[uint i] { get { if (i >= this.sa.Length) lock (this.lck) { var lngth = this.sa.Length; while (i >= lngth) { var bf = (ushort[])MCPY.Clone(); if (lngth == 0) { for (uint bi = 0, wi = 0, w = 0, msk = 0x8000, v = 0; w < bf.Length; bi += WHLPTRN[wi++], wi = (wi >= WHTS) ? 0 : wi) { if (msk >= 0x8000) { msk = 1; v = bf[w++]; } else msk <<= 1; if ((v & msk) == 0) { var p = FSTBP + (bi + bi); var k = (p * p - FSTBP) >> 1; if (k >= PGRNG) break; var pd = p / WCRC; var kd = k / WCRC; var kn = WHLNDX[k - kd * WCRC]; for (uint wrd = kd * WPC + (uint)(kn >> 4), ndx = wi * WHTS + kn; wrd < bf.Length; ) { var st = WSLUT[ndx]; bf[wrd] |= st.msk; wrd += st.mlt * pd + st.xtr; ndx = st.nxt; } } } } else { this.lwi += PGRNG; cullbf(this.lwi, bf); } var c = count(PGRNG, bf); var na = new byte[lngth + c]; sa.CopyTo(na, 0); for (uint p = FSTBP + (this.lwi << 1), wi = 0, w = 0, msk = 0x8000, v = 0; lngth < na.Length; p += (uint)(WHLPTRN[wi++] << 1), wi = (wi >= WHTS) ? 0 : wi) { if (msk >= 0x8000) { msk = 1; v = bf[w++]; } else msk <<= 1; if ((v & msk) == 0) { var pd = p / WCRC; na[lngth++] = (byte)(((pd - this.lpd) << 6) + wi); this.lpd = pd; } } this.sa = na; } } return this.sa[i]; } } } static readonly Bpa baseprms = new Bpa(); //the base primes array using the above class struct PrcsSpc { public Task tsk; public ushort[] buf; } //used for multi-threading buffer array processing #endregion #region private static methods static int count(uint bitlim, ushort[] buf) { //very fast counting method using the CLUT look up table if (bitlim < BFRNG) { var addr = (bitlim - 1) / WCRC; var bit = WHLNDX[bitlim - addr * WCRC] - 1; addr *= WPC; for (var i = 0; i < 3; ++i) buf[addr++] |= (ushort)((unchecked((ulong)-2) << bit) >> (i << 4)); } var acc = 0; for (uint i = 0, w = 0; i < bitlim; i += WCRC) acc += CLUT[buf[w++]] + CLUT[buf[w++]] + CLUT[buf[w++]]; return acc; } static void cullbf(ulong lwi, ushort[] b) { //fast buffer segment culling method using a Wheel State Look Up Table ulong nlwi = lwi; for (var i = 0u; i < b.Length; nlwi += PGRNG, i += PGSZ) MCPY.CopyTo(b, i); //copy preculled lower base primes. for (uint i = 0, pd = 0; ; ++i) { pd += (uint)baseprms[i] >> 6; var wi = baseprms[i] & 0x3Fu; var wp = (uint)WHLPOS[wi]; var p = pd * WCRC + PRLUT[wi]; var k = ((ulong)p * (ulong)p - FSTBP) >> 1; if (k >= nlwi) break; if (k < lwi) { k = (lwi - k) % (WCRC * p); if (k != 0) { var nwp = wp + (uint)((k + p - 1) / p); k = (WHLRNDUP[nwp] - wp) * p - k; } } else k -= lwi; var kd = k / WCRC; var kn = WHLNDX[k - kd * WCRC]; for (uint wrd = (uint)kd * WPC + (uint)(kn >> 4), ndx = wi * WHTS + kn; wrd < b.Length; ) { var st = WSLUT[ndx]; b[wrd] |= st.msk; wrd += st.mlt * pd + st.xtr; ndx = st.nxt; } } } static Task cullbftsk(ulong lwi, ushort[] b, Action f) { // forms a task of the cull buffer operaion return Task.Factory.StartNew(() => { cullbf(lwi, b); f(b); }); } //iterates the action over each page up to the page including the top_number, //making an adjustment to the top limit for the last page. //this method works for non-dependent actions that can be executed in any order. static void IterateTo(ulong top_number, Action actn) { PrcsSpc[] ps = new PrcsSpc[NUMPRCSPCS]; for (var s = 0u; s < NUMPRCSPCS; ++s) ps[s] = new PrcsSpc { buf = new ushort[BFSZ], tsk = Task.Factory.StartNew(() => { }) }; var topndx = (top_number - FSTBP) >> 1; for (ulong ndx = 0; ndx <= topndx; ) { ps[0].tsk.Wait(); var buf = ps[0].buf; for (var s = 0u; s < NUMPRCSPCS - 1; ++s) ps[s] = ps[s + 1]; var lowi = ndx; var nxtndx = ndx + BFRNG; var lim = topndx < nxtndx ? (uint)(topndx - ndx + 1) : BFRNG; ps[NUMPRCSPCS - 1] = new PrcsSpc { buf = buf, tsk = cullbftsk(ndx, buf, (b) => actn(lowi, lim, b)) }; ndx = nxtndx; } for (var s = 0u; s < NUMPRCSPCS; ++s) ps[s].tsk.Wait(); } //iterates the predicate over each page up to the page where the predicate paramenter returns true, //this method works for dependent operations that need to be executed in increasing order. //it is somewhat slower than the above as the predicate function is executed outside the task. static void IterateUntil(Func prdct) { PrcsSpc[] ps = new PrcsSpc[NUMPRCSPCS]; for (var s = 0u; s < NUMPRCSPCS; ++s) { var buf = new ushort[BFSZ]; ps[s] = new PrcsSpc { buf = buf, tsk = cullbftsk(s * BFRNG, buf, (bfr) => { }) }; } for (var ndx = 0UL; ; ndx += BFRNG) { ps[0].tsk.Wait(); var buf = ps[0].buf; var lowi = ndx; if (prdct(lowi, buf)) break; for (var s = 0u; s < NUMPRCSPCS - 1; ++s) ps[s] = ps[s + 1]; ps[NUMPRCSPCS - 1] = new PrcsSpc { buf = buf, tsk = cullbftsk(ndx + NUMPRCSPCS * BFRNG, buf, (bfr) => { }) }; } } #endregion #region initialization /// /// the static constructor is used to initialize the static readonly arrays. /// static UltimatePrimesSoE() { WHLPOS = new byte[WHLPTRN.Length + 1]; //to look up wheel position index from wheel index for (byte i = 0, acc = 0; i < WHLPTRN.Length; ++i) { acc += WHLPTRN[i]; WHLPOS[i + 1] = acc; } WHLNDX = new byte[WCRC + 1]; for (byte i = 1; i < WHLPOS.Length; ++i) { for (byte j = (byte)(WHLPOS[i - 1] + 1); j <= WHLPOS[i]; ++j) WHLNDX[j] = i; } WHLRNDUP = new byte[WCRC * 2]; for (byte i = 1; i < WHLRNDUP.Length; ++i) { if (i > WCRC) WHLRNDUP[i] = (byte)(WCRC + WHLPOS[WHLNDX[i - WCRC]]); else WHLRNDUP[i] = WHLPOS[WHLNDX[i]]; } Func nmbts = (v) => { var acc = 0; while (v != 0) { acc += (int)v & 1; v >>= 1; } return acc; }; CLUT = new byte[1 << 16]; for (var i = 0; i < CLUT.Length; ++i) CLUT[i] = (byte)nmbts((ushort)(i ^ -1)); PRLUT = new byte[WHTS]; for (var i = 0; i < PRLUT.Length; ++i) { var t = (uint)(WHLPOS[i] * 2) + FSTBP; if (t >= WCRC) t -= WCRC; if (t >= WCRC) t -= WCRC; PRLUT[i] = (byte)t; } WSLUT = new Wst[WHTS * WHTS]; for (var x = 0u; x < WHTS; ++x) { var p = FSTBP + 2u * WHLPOS[x]; var pr = p % WCRC; for (uint y = 0, pos = (p * p - FSTBP) / 2; y < WHTS; ++y) { var m = WHLPTRN[(x + y) % WHTS]; pos %= WCRC; var posn = WHLNDX[pos]; pos += m * pr; var nposd = pos / WCRC; var nposn = WHLNDX[pos - nposd * WCRC]; WSLUT[x * WHTS + posn] = new Wst { msk = (ushort)(1 << (int)(posn & 0xF)), mlt = (byte)(m * WPC), xtr = (byte)(WPC * nposd + (nposn >> 4) - (posn >> 4)), nxt = (ushort)(WHTS * x + nposn) }; } } MCPY = new ushort[PGSZ]; foreach (var lp in BWHLPRMS.SkipWhile(p => p < FSTCP)) { var p = (uint)lp; var k = (p * p - FSTBP) >> 1; var pd = p / WCRC; var kd = k / WCRC; var kn = WHLNDX[k - kd * WCRC]; for (uint w = kd * WPC + (uint)(kn >> 4), ndx = WHLNDX[(2 * WCRC + p - FSTBP) / 2] * WHTS + kn; w < MCPY.Length; ) { var st = WSLUT[ndx]; MCPY[w] |= st.msk; w += st.mlt * pd + st.xtr; ndx = st.nxt; } } } #endregion #region public class // this class implements the enumeration (IEnumerator). // it works by farming out tasks culling pages, which it then processes in order by // enumerating the found primes as recognized by the remaining non-composite bits // in the cull page buffers. class nmrtr : IEnumerator , IEnumerator, IDisposable { PrcsSpc[] ps = new PrcsSpc[NUMPRCSPCS]; ushort[] buf; public nmrtr() { for (var s = 0u; s < NUMPRCSPCS; ++s) ps[s] = new PrcsSpc { buf = new ushort[BFSZ] }; for (var s = 1u; s < NUMPRCSPCS; ++s) { ps[s].tsk = cullbftsk((s - 1u) * BFRNG, ps[s].buf, (bfr) => { }); } buf = ps[0].buf; } ulong _curr, i = (ulong)-WHLPTRN[WHTS - 1]; int b = -BWHLPRMS.Length - 1; uint wi = WHTS - 1; ushort v, msk = 0; public ulong Current { get { return this._curr; } } object IEnumerator.Current { get { return this._curr; } } public bool MoveNext() { if (b < 0) { if (b == -1) b += buf.Length; //no yield!!! so automatically comes around again else { this._curr = (ulong)BWHLPRMS[BWHLPRMS.Length + (++b)]; return true; } } do { i += WHLPTRN[wi++]; if (wi >= WHTS) wi = 0; if ((this.msk <<= 1) == 0) { if (++b >= BFSZ) { b = 0; for (var prc = 0; prc < NUMPRCSPCS - 1; ++prc) ps[prc] = ps[prc + 1]; ps[NUMPRCSPCS - 1u].buf = buf; ps[NUMPRCSPCS - 1u].tsk = cullbftsk(i + (NUMPRCSPCS - 1u) * BFRNG, buf, (bfr) => { }); ps[0].tsk.Wait(); buf = ps[0].buf; } v = buf[b]; this.msk = 1; } } while ((v & msk) != 0u); _curr = FSTBP + i + i; return true; } public void Reset() { throw new Exception("Primes enumeration reset not implemented!!!"); } public void Dispose() { } } #endregion #region public instance method and associated sub private method /// /// Gets the enumerator for the primes. /// /// The enumerator of the primes. public IEnumerator GetEnumerator() { return new nmrtr(); } /// /// Gets the enumerator for the primes. /// /// The enumerator of the primes. IEnumerator IEnumerable.GetEnumerator() { return new nmrtr(); } #endregion #region public static methods /// /// Gets the count of primes up the number, inclusively. /// /// The ulong top number to check for prime. /// The long number of primes found. public static long CountTo(ulong top_number) { if (top_number < FSTBP) return BWHLPRMS.TakeWhile(p => p <= top_number).Count(); var cnt = (long)BWHLPRMS.Length; IterateTo(top_number, (lowi, lim, b) => { Interlocked.Add(ref cnt, count(lim, b)); }); return cnt; } /// /// Gets the sum of the primes up the number, inclusively. /// /// The uint top number to check for prime. /// The ulong sum of all the primes found. public static ulong SumTo(uint top_number) { if (top_number < FSTBP) return (ulong)BWHLPRMS.TakeWhile(p => p <= top_number).Aggregate(0u, (acc, p) => acc += p); var sum = (long)BWHLPRMS.Aggregate(0u, (acc, p) => acc += p); Func sumbf = (lowi, bitlim, buf) => { var acc = 0L; for (uint i = 0, wi = 0, msk = 0x8000, w = 0, v = 0; i < bitlim; i += WHLPTRN[wi++], wi = wi >= WHTS ? 0 : wi) { if (msk >= 0x8000) { msk = 1; v = buf[w++]; } else msk <<= 1; if ((v & msk) == 0) acc += (long)(FSTBP + ((lowi + i) << 1)); } return acc; }; IterateTo(top_number, (pos, lim, b) => { Interlocked.Add(ref sum, sumbf(pos, lim, b)); }); return (ulong)sum; } /// /// Gets the prime number at the zero based index number given. /// /// The long zero-based index number for the prime. /// The ulong prime found at the given index. public static ulong ElementAt(long index) { if (index < BWHLPRMS.Length) return (ulong)BWHLPRMS.ElementAt((int)index); long cnt = BWHLPRMS.Length; var ndx = 0UL; var cycl = 0u; var bit = 0u; IterateUntil((lwi, bfr) => { var c = count(BFRNG, bfr); if ((cnt += c) < index) return false; ndx = lwi; cnt -= c; c = 0; do { var w = cycl++ * WPC; c = CLUT[bfr[w++]] + CLUT[bfr[w++]] + CLUT[bfr[w]]; cnt += c; } while (cnt < index); cnt -= c; var y = (--cycl) * WPC; ulong v = ((ulong)bfr[y + 2] << 32) + ((ulong)bfr[y + 1] << 16) + bfr[y]; do { if ((v & (1UL << ((int)bit++))) == 0) ++cnt; } while (cnt <= index); --bit; return true; }); return FSTBP + ((ndx + cycl * WCRC + WHLPOS[bit]) << 1); } #endregion }
上面的代码将在四核（8个线程，包括HT）i7-2700K（3.5 GHz）上计算大约1.55秒的质数到10亿次，而你的E7500的速度可能会慢4倍，因为线程更少，时钟更少速度。 大约四分之三的时间是运行枚举MoveNext（）方法和Current属性的时间，所以我提供公共静态方法“CountTo”，“SumTo”和“ElementAt”来计算一个素数的数量或总和。范围和第n个基于零的素数，分别不使用枚举。 使用UltimatePrimesSoE.CountTo（1000000000）静态方法在我的机器上在大约0.32秒内产生50847534，因此在英特尔E7500上的时间不应超过大约1.28秒。
EDIT_ADD：有趣的是，这个代码在x86 32位模式下的运行速度比在x64 64位模式下快30％，这可能是因为避免了将uint32数字扩展到ulong的额外开销。 以上所有时序均适用于64位模式。 END_EDIT_ADD
在近300行（密集）代码行中，这种实现并不简单，但这是执行所有所述优化的成本，这些优化使得此代码非常有效。 Aaron Murgatroyd的另一个答案并不是那么多的代码行; 虽然他的代码密度较低，但他的代码也慢了四倍。 实际上，几乎所有的执行时间都花费在我的代码的私有静态“cullbf”方法的最后“for循环”中，该方法只有四个语句长加上范围条件检查; 所有其余的代码只是支持该循环的重复应用。
这个代码比其他答案更快的原因是这个代码比你的代码更快的原因，而不是他只进行奇数素数候选者的步骤（1）优化。 他使用多处理几乎完全没有效率，因为只有30％的优势，而不是真正的四核CPU应该可能的四倍因子，因为它是每个素数而不是小页面上的所有素数，并且使用不安全的指针数组访问作为消除每个循环的数组绑定检查的DotNet计算成本的方法实际上减慢了代码与直接使用包含边界检查的数组相比，因为DotNet Just In Time（JIT）编译器产生非常低效的代码用于指针访问。 另外，他的代码枚举了素数就像我的代码可以做的那样，枚举的每个枚举素数有10个CPU时钟周期成本，在他的情况下也略差一些，因为他使用了内置的C＃迭代器，比我的“自己动手”的IEnumerator界面更有效率。 但是，为了获得最大速度，我们应该完全避免枚举; 然而，即使他提供的“Count”实例方法也使用“foreach”循环，这意味着枚举。
总而言之，这个答案代码产生的主要答案比E7500 CPU上的问题代码快25倍（在具有更多内核/线程的CPU上快多倍）使用更少的内存，并且不限于更小的主要范围。 32位数字范围，但代价是代码复杂性增加。
我的multithreading实现（需要.NET 4.0）：
using System; using System.Collections; using System.Collections.Generic; using System.Threading.Tasks; namespace PrimeGenerator { // The block element type for the bit array, // use any unsigned value. WARNING: UInt64 is // slower even on x64 architectures. using BitArrayType = System.UInt32; // This should never be any bigger than 256 bits - leave as is. using BitsPerBlockType = System.Byte; // The prime data type, this can be any unsigned value, the limit // of this type determines the limit of Prime value that can be // found. WARNING: UInt64 is slower even on x64 architectures. using PrimeType = System.UInt32; /// /// Calculates prime number using the Sieve of Eratosthenes method. /// /// /// /// var lpPrimes = new Eratosthenes(1e7); /// foreach (UInt32 luiPrime in lpPrimes) /// Console.WriteLine(luiPrime); /// public class Eratosthenes : IEnumerable { #region Constants /// /// Constant for number of bits per block, calculated based on size of BitArrayType. /// const BitsPerBlockType cbBitsPerBlock = sizeof(BitArrayType) * 8; #endregion #region Protected Locals /// /// The limit for the maximum prime value to find. /// protected readonly PrimeType mpLimit; /// /// True if the class is multi-threaded /// protected readonly bool mbThreaded; /// /// The current bit array where a set bit means /// the odd value at that location has been determined /// to not be prime. /// protected BitArrayType[] mbaOddNotPrime; #endregion #region Initialisation /// /// Create Sieve of Eratosthenes generator. /// /// The limit for the maximum prime value to find. /// True if threaded, false otherwise. public Eratosthenes(PrimeType limit, bool threaded) { // Check limit range if (limit > PrimeType.MaxValue - (PrimeType)Math.Sqrt(PrimeType.MaxValue)) throw new ArgumentOutOfRangeException(); mbThreaded = threaded; mpLimit = limit; FindPrimes(); } /// /// Create Sieve of Eratosthenes generator in multi-threaded mode. /// /// The limit for the maximum prime value to find. public Eratosthenes(PrimeType limit) : this(limit, true) { } #endregion #region Private Methods /// /// Calculates compartment indexes for a multi-threaded operation. /// /// The inclusive starting index. /// The exclusive ending index. /// The number of threads. /// An array of thread elements plus 1 containing the starting and exclusive ending indexes to process for each thread. private PrimeType[] CalculateCompartments(PrimeType startInclusive, PrimeType endExclusive, ref int threads) { if (threads == 0) threads = 1; if (threads == -1) threads = Environment.ProcessorCount; if (threads > endExclusive - startInclusive) threads = (int)(endExclusive - startInclusive); PrimeType[] liThreadIndexes = new PrimeType[threads + 1]; liThreadIndexes[threads] = endExclusive; PrimeType liIndexesPerThread = (endExclusive - startInclusive) / (PrimeType)threads; for (PrimeType liCount = 0; liCount < threads; liCount++) liThreadIndexes[liCount] = liCount * liIndexesPerThread; return liThreadIndexes; } /// /// Executes a simple for loop in parallel but only creates /// a set amount of threads breaking the work up evenly per thread, /// calling the body only once per thread, this is different /// to the .NET 4.0 For method which calls the body for each index. /// /// The type of parameter to pass to the body. /// The starting index. /// The exclusive ending index. /// The parameter to pass to the body. /// The body to execute per thread. /// The number of threads to execute. private void For( PrimeType startInclusive, PrimeType endExclusive, ParamType parameter, Action body, int threads) { PrimeType[] liThreadIndexes = CalculateCompartments(startInclusive, endExclusive, ref threads); if (threads > 1) Parallel.For( 0, threads, new System.Threading.Tasks.ParallelOptions(), (liThread) => { body(liThreadIndexes[liThread], liThreadIndexes[liThread + 1], parameter); } ); else body(startInclusive, endExclusive, parameter); } /// /// Finds the prime number within range. /// private unsafe void FindPrimes() { // Allocate bit array. mbaOddNotPrime = new BitArrayType[(((mpLimit >> 1) + 1) / cbBitsPerBlock) + 1]; // Cache Sqrt of limit. PrimeType lpSQRT = (PrimeType)Math.Sqrt(mpLimit); int liThreads = Environment.ProcessorCount; if (!Threaded) liThreads = 0; // Fix the bit array for pointer access fixed (BitArrayType* lpbOddNotPrime = &mbaOddNotPrime[0]) { IntPtr lipBits = (IntPtr)lpbOddNotPrime; // Scan primes up to lpSQRT for (PrimeType lpN = 3; lpN <= lpSQRT; lpN += 2) { // If the current bit value for index lpN is cleared (prime) if ( ( lpbOddNotPrime[(lpN >> 1) / cbBitsPerBlock] & ((BitArrayType)1 << (BitsPerBlockType)((lpN >> 1) % cbBitsPerBlock)) ) == 0 ) { // If multi-threaded if (liThreads > 1) { // Leave it cleared (prime) and mark all multiples of lpN*2 from lpN*lpN as not prime For( 0, ((mpLimit - (lpN * lpN)) / (lpN << 1)) + 1, lpN, (start, end, n) => { BitArrayType* lpbBits = (BitArrayType*)lipBits; PrimeType lpM = n * n + (start * (n << 1)); for (PrimeType lpCount = start; lpCount < end; lpCount++) { // Set as not prime lpbBits[(lpM >> 1) / cbBitsPerBlock] |= (BitArrayType)((BitArrayType)1 << (BitsPerBlockType)((lpM >> 1) % cbBitsPerBlock)); lpM += n << 1; } }, liThreads); } else { // Leave it cleared (prime) and mark all multiples of lpN*2 from lpN*lpN as not prime for (PrimeType lpM = lpN * lpN; lpM <= mpLimit; lpM += lpN<<1) // Set as not prime lpbOddNotPrime[(lpM >> 1) / cbBitsPerBlock] |= (BitArrayType)((BitArrayType)1 << (BitsPerBlockType)((lpM >> 1) % cbBitsPerBlock)); } } } } } /// /// Gets a bit value by index. /// /// The blocks containing the bits. /// The index of the bit. /// True if bit is set, false if cleared. private bool GetBitSafe(BitArrayType[] bits, PrimeType index) { return (bits[index / cbBitsPerBlock] & ((BitArrayType)1 << (BitsPerBlockType)(index % cbBitsPerBlock))) != 0; } #endregion #region Public Properties /// /// Gets whether this class is multi-threaded or not. /// public bool Threaded { get { return mbThreaded; } } /// /// Get the limit for the maximum prime value to find. /// public PrimeType Limit { get { return mpLimit; } } /// /// Returns the number of primes found in the range. /// public PrimeType Count { get { PrimeType lptCount = 0; foreach (PrimeType liPrime in this) lptCount++; return lptCount; } } /// /// Determines if a value in range is prime or not. /// /// The value to test for primality. /// True if the value is prime, false otherwise. public bool this[PrimeType test] { get { if (test > mpLimit) throw new ArgumentOutOfRangeException(); if (test <= 1) return false; if (test == 2) return true; if ((test & 1) == 0) return false; return !GetBitSafe(mbaOddNotPrime, test >> 1); } } #endregion #region Public Methods /// /// Gets the enumerator for the primes. /// /// The enumerator of the primes. public IEnumerator GetEnumerator() { // Two always prime. yield return 2; // Start at first block, second MSB. int liBlock = 0; byte lbBit = 1; BitArrayType lbaCurrent = mbaOddNotPrime[0] >> 1; // For each value in range stepping in incrments of two for odd values. for (PrimeType lpN = 3; lpN <= mpLimit; lpN += 2) { // If current bit not set then value is prime. if ((lbaCurrent & 1) == 0) yield return lpN; // Move to NSB. lbaCurrent >>= 1; // Increment bit value. lbBit++; // If block is finished. if (lbBit == cbBitsPerBlock) { // Move to first bit of next block. lbBit = 0; liBlock++; lbaCurrent = mbaOddNotPrime[liBlock]; } } } #endregion #region IEnumerable Implementation /// /// Gets the enumerator for the primes. /// ///
The multi-threading works by threading the inner most loop, this way there are no data locking issues because the multiple threads work with a subset of the array and dont overlap for each job done.
Seems to be quite fast, can generate all primes up to a limit of 1,000,000,000 on an AMD Phenom II X4 965 processor in 5.8 seconds. Special implementations like the Atkins are faster, but this is fast for the Sieve of Eratosthenes.
A while back I tried to implement The Sieve of Atkin in parallell . It was a failure. I haven’t done any deeper research but it seems that both Sieve Of Eratosthenes and The Sieve of Atkin are hard to scale over multiple CPUs because the implementations I’ve seen uses a list of integers that is shared. Shared state is a heavy anchor to carry when you try to scale over multiple CPUs.