For primality testing of 64 bit ulong, I have optimized a very fast trial-by-division test using possible factors of the form 6k+/-1.  For input numbers less than uint.MaxValue serial testing is performed.  For really larger numbers, parallel processing is performed.Again as optimized as it is, it is still trial-by-division and does not employ any Miller Rabin techniques.I have 2 signature overloads: one is for long which works with all signed value type integers, and the other is for ulong which works with the unsigned value type integers, i.e. everything but BigInteger since it is not a value type.public static bool IsPrime(long number){    if (number < 2) { return false; }    return IsPrime((ulong)number);}public static bool IsPrime(ulong number){    // Get the quick checks out of the way.    if (number < 2) { return false; }    // Dispense with multiples of 2 and 3.    if (number % 2 == 0) { return (number == 2); }    if (number % 3 == 0) { return (number == 3); }    // Another quick check to eliminate known composites.    // http://programmers.stackexchange.com/questions/120934/best-and-most-used-algorithm-for-finding-the-primality-of-given-positive-number/120963#120963    if (!(((number - 1) % 6 == 0) || ((number + 1) % 6 == 0)))    {        return false;    }    // Quick checks are over.  What remains is a POSSIBLE prime.    // Must iterate to determine absolute the answer.    // We loop over 1/6 of the required possible factors to check,    // but since we check twice in each iteration, we are actually    // checking 1/3 of the possible divisors.  This is an improvement    // over the typical naive test of odds only which tests 1/2    // of the factors.    // Though the whole number portion of the square root of ulong.MaxValue    // would fit in a uint, we want to cast to uint and then back to ulong.    ulong root = (ulong)(uint)Math.Sqrt(number);    // Fix Corner Case: Math.Sqrt error for really HUGE ulong.    if (root == 0) root = (ulong)uint.MaxValue;    // For small enough numbers, serial is faster than parallel.    // Obviously there is some number where parallel becomes faster.    // I do not know at which point that occurs.    // I have arbitrarily chosen a point based on the square root.    //     rootCutoff = 33,554,432    //     square of rootCutoff = 1,125,899,906,842,624    const ulong rootCutoff = 65536UL * 512UL;    if ((root < rootCutoff) || (Environment.ProcessorCount == 1))    {        // Serial Loop for smaller numbers:        // Start at 5, which is (6k-1) where k=1.        // Increment the loop by 6, which is same as incrementing k by 1.        for (ulong factor = 5; factor <= root; factor += 6)        {            // Check (6k-1)            if (number % factor == 0) { return false; }            // Check (6k+1)            if (number % (factor + 2) == 0) { return false; }        }        return true;    }    // Parallel Looping for the bigger numbers:    return IsPrimeParallel(number, root);}private static bool IsPrimeParallel(ulong number, ulong root){    int composite = 0;    // I arbitrarily choose the number of chunks to be the same as the number of processors.    // Each chunk will processed in its own thread, but this in no way is equivalent    // of saying each thread goes to its own core, or vice versa.    int chunks = Environment.ProcessorCount;    // perform cast once rather than a billion times    ulong chunks64 = (ulong)chunks;    int kEnd = (int)(root / 6) + 1;    int iEnd = (kEnd / chunks) + 1;    Parallel.For(0, chunks, (chunk, loopState) =>    {        // perform cast once rather than 90-715 million times        ulong offset = (ulong)chunk + 1;        for (int i = 0; i < iEnd; i++)        {            ulong k = (chunks64 * (ulong)i) + offset;            if (k > root) { break; }            ulong factor = (6 * k) - 1;            if (number % factor == 0) // (6k-1)            {                Interlocked.Exchange(ref composite, 1);                loopState.Stop();            }            else if (number % (factor + 2) == 0) // (6k+1)            {                Interlocked.Exchange(ref composite, 1);                loopState.Stop();            }            if (loopState.IsStopped) { break; }        }    });    return (composite == 0);}Avoid Too Many ThreadsAs Ive seen on SO and CR, parallel processing can easily be done wrong and take much longer than simple serial.  Inside my loop, there are only a few, fast calculations.  It would be a performance drag to create a thread for each k (worst case is over 700 million).This normally is a great candidate for a range Partitioner but I need more than a simple range.Consider a simplified example where the domain of k is 1 to 100 inclusively, which is to say that kEnd is 101.  I also want to generate no more than 10 threads in this example.I do NOT want ranges like:Range 0 is {1, 2, 3 , , 8, 9, 10}Range 1 is {11, 12, 13 , , 18,    19, 20} . . . Range 9 is {91, 92, 93 , , 98, 99, 100}I want series like:Series 0 is {1, 11, 21, , 71, 81, 91}Series 1 is {2, 12, 22, , 72, 82, 92}. . . Series 9 is {10, 20, 30, , 80, 90, 100}Square Root Corner CaseThere is a special corner case where Math.Sqrt(ulong.MaxValue) returns the wrong value.  This is not the fault of Sqrt itself but just the nature of the (implicit) cast of ulong.MaxValue to a double since you have an integer value fully and exactly represented by 64 bits but you are squeezing them into a 64 bit floating point approximation.This corner case is easy to detect and just as easy to correct.  The key is to be aware of it in the first place.Performance[Edit: I erroneously listed the 31 & 32 bit times as seconds.  The correct unit is milliseconds.]On my 8-core laptop in serial only mode:Largest 31 bit prime takes 0.17 milliseconds.Largest 32 bit prime takes 0.23 milliseconds.Largest 63 bit prime takes over 13 seconds.Largest 64 bit prime takes over 18 seconds.In parallel mode:Largest 63 bit prime takes 6.6 seconds.Largest 64 bit prime takes 9.28 seconds.QuestionsBeing this is CR, there is always an implied question of Do you have any constructive comments?While I have used Parallel.For many times before, this is my first implementation where I had to create a specifically arranged series rather than a simple range.  Is this done correctly and/or could it be done better?Other than using Miller-Rabin techniques, can this be made faster?  Ive used more threads and less threads - or really what I call chunks and assuming that each chunk gets its own thread - but on my laptop the fastest times consistently were when I used a chunk count equal to my processor count.
Optimized ulong prime test using 6k+/-1 in parallel threads with C#
c#;primes;task parallel library
If I understand correctly, the (k > root) test may only trigger on a few values at the end of the loop to avoid computation.  If it is indeed the case then it may be cheaper to avoid the test altogether and risk doing more work.--I would try to have a loop on factor and throw i and k away : no multiplications, less additions.I quickly did a test that approximates the original loop to have an idea of the gain:    public static ulong IsPrimeInternalLoop(ulong number, ulong root)    {        ulong ret = 0;        int composite = 0;        // I arbitrarily choose the number of chunks to be the same as the number of processors.        // Each chunk will processed in its own thread, but this in no way is equivalent        // of saying each thread goes to its own core, or vice versa.        int chunks = Environment.ProcessorCount;        // perform cast once rather than a billion times        ulong chunks64 = (ulong)chunks;        int kEnd = (int)(root / 6) + 1;        int iEnd = (kEnd / chunks) + 1;        for(int chunk = 0; chunk < chunks; ++chunk)        {            // perform cast once rather than 90-715 million times            ulong offset = (ulong)chunk + 1;            for(int i = 0; i < iEnd; i++)            {                ulong k = (chunks64 * (ulong)i) + offset;                if(k > root) { break; }                ulong factor = (6 * k) - 1;                ++ret;            }        }        return ret;    }// without i:    public static ulong IsPrimeInternalLoop2(ulong number, ulong root)    {        ulong ret = 0;        int composite = 0;        // I arbitrarily choose the number of chunks to be the same as the number of processors.        // Each chunk will processed in its own thread, but this in no way is equivalent        // of saying each thread goes to its own core, or vice versa.        int chunks = Environment.ProcessorCount;        // perform cast once rather than a billion times        ulong chunks64 = (ulong)chunks;        //int kEnd = (int)(root / 6) + 1;        ulong kEnd = (root / 6) + chunks64;        for(int chunk = 0; chunk < chunks; ++chunk)        {            for(ulong k = (ulong)chunk + 1; k <= kEnd; k += chunks64)            {                ulong factor = (6 * k) - 1;                ++ret;            }        }        return ret;    }// without k:    public static ulong IsPrimeInternalLoop3(ulong number, ulong root)    {        ulong ret = 0;        int composite = 0;        // I arbitrarily choose the number of chunks to be the same as the number of processors.        // Each chunk will processed in its own thread, but this in no way is equivalent        // of saying each thread goes to its own core, or vice versa.        int chunks = Environment.ProcessorCount;        // perform cast once rather than a billion times        ulong chunks64 = (ulong)chunks;        ulong kEnd = (root / 6) + chunks64;        for(int chunk = 0; chunk < chunks; ++chunk)        {            ulong fEnd = (6 * kEnd) - 1;            ulong fStep = (6 * chunks64);            // perform cast once rather than 90-715 million times            for(ulong factor = 6 * ((ulong)chunk + 1) - 1; factor <= fEnd; factor += fStep)            {                ++ret;            }        }        return ret;    }Note: functions 2 & 3 do a bit too many tests.  It would need to be more accurate. Still, the idea is there.I tested the above with the biggest prime on 64 bits and its neighbors:        ulong[] values = new ulong[] { 18446744073709551557 - 1, 18446744073709551557, 18446744073709551557 + 1 };        foreach(ulong v in values)        {            ulong root = (ulong)Math.Sqrt(v);            if(root <= 0)            {                root = int.MaxValue;            }            long d0 = long.MaxValue, d1 = long.MaxValue , d2 = long.MaxValue;            ulong a = 0, b = 0, c = 0;            for(int i = 0; i < 5; i++)            {                long start = DateTime.Now.Ticks;                a = Form1.IsPrimeInternalLoop(v, root);                long stop = DateTime.Now.Ticks;                d0 = Math.Min(stop - start, d0);            }            for(int i = 0; i < 5; i++)            {                long start = DateTime.Now.Ticks;                b = Form1.IsPrimeInternalLoop2(v, root);                long stop = DateTime.Now.Ticks;                d1 = Math.Min(stop - start, d1);            }            for(int i = 0; i < 5; i++)            {                long start = DateTime.Now.Ticks;                c = Form1.IsPrimeInternalLoop3(v, root);                long stop = DateTime.Now.Ticks;                d2 = Math.Min(stop - start, d2);            }            Console.WriteLine(value: {0} root: {1}, v, root);            Console.WriteLine(times: {0} {1} {2}, d0, d1, d2);            Console.WriteLine(counts: {0} {1} {2}, a, b, c);        }The output shows smaller times on both modifications (broadly 9%) (times)The number of factors computed shows slightly more work on functions 2 and 3 than on the original one (counts)value: 18446744073709551556 root: 4294967296times: 15132027 13884024 13728024counts: 715827888 715827889 715827890value: 18446744073709551557 root: 4294967296times: 15288027 13884024 13728024counts: 715827888 715827889 715827890value: 18446744073709551558 root: 4294967296times: 15132027 13884024 13728024counts: 715827888 715827889 715827890--But this is only loop cost, and against the multithread cost, the above is nothing.Calling loopState.IsStopped may be a bit costly compared to the computations.  So you may want to reduce it's usage, for example: having an arbitrary countdown only testing for it every few iterations.This time I ran IsPrime in a loop. (i7 2600K/W7 VS15 release/x64)The original code time: ~31000000With local countdown of 4k around IsStopped: ~25000000 (I actually underestimated the gain)When combining the above with the loop of the function #3: ~21000000So here is the fastest code I got.  (You still need to fix the loop count)    private static bool IsPrimeParallel(ulong number, ulong root)    {        int composite = 0;        // I arbitrarily choose the number of chunks to be the same as the number of processors.        // Each chunk will processed in its own thread, but this in no way is equivalent        // of saying each thread goes to its own core, or vice versa.        int chunks = Environment.ProcessorCount;        // perform cast once rather than a billion times        Parallel.For(0, chunks, (chunk, loopState) =>        {            // perform cast once rather than 90-715 million times            int iscd = 4096; // is stopping count down            // access of variables outside of the lambda function may be slow, depending on how the compiler does it            ulong chunks64 = (ulong)chunks;            ulong kEnd = (root / 6) + chunks64;            ulong fEnd = (6 * kEnd) - 1;            ulong fStep = (6 * chunks64);            // perform cast once rather than 90-715 million times            for(ulong factor = 6 * ((ulong)chunk + 1) - 1; factor <= fEnd; factor += fStep)            {                if(number % factor == 0) // (6k-1)                {                    Interlocked.Exchange(ref composite, 1);                    loopState.Stop();                }                else if(number % (factor + 2) == 0) // (6k+1)                {                    Interlocked.Exchange(ref composite, 1);                    loopState.Stop();                }                if(--iscd == 0)                {                    if(loopState.IsStopped) { break; }                    iscd = 4096;                }            }        });        return (composite == 0);    }