(*^ ::[ frontEndVersion = "Microsoft Windows Mathematica Notebook Front End Version 2.2"; microsoftWindowsStandardFontEncoding; fontset = title, "Times", 24, L0, center, nohscroll, bold; fontset = subtitle, "Times", 18, L0, center, nohscroll, bold; fontset = subsubtitle, "Times", 14, L0, center, nohscroll, italic; fontset = section, "Times", 18, L0, nohscroll, bold, grayBox; fontset = subsection, "Times", 14, L0, nohscroll, bold, blackBox; fontset = subsubsection, "Times", 12, L0, nohscroll, bold, whiteBox; fontset = text, "Times", 12, L0, nohscroll; fontset = smalltext, "Times", 10, L0, nohscroll; fontset = input, "Courier", 12, L-5, nowordwrap, bold; fontset = output, "Courier", 12, L-5, nowordwrap; fontset = message, "Courier", 12, L-5, nowordwrap, R65535; fontset = print, "Courier", 12, L-5, nowordwrap; fontset = info, "Courier", 12, L-5, nowordwrap, B65535; fontset = postscript, "Courier", 12, L0, nowordwrap; fontset = name, "Geneva", 10, L0, nohscroll, italic; fontset = header, "Times", 12, L0; fontset = footer, "Times", 12, L0, center; fontset = help, "Times", 10, L0, nohscroll; fontset = clipboard, "Times", 12, L0, nohscroll; fontset = completions, "Times", 12, L0, nohscroll; fontset = graphics, "Courier New", 10, L0, nowordwrap, nohscroll; fontset = special1, "Times", 12, L0, nohscroll; fontset = special2, "Times", 12, L0, nohscroll; fontset = special3, "Times", 12, L0, nohscroll; fontset = special4, "Times", 12, L0, nohscroll; fontset = special5, "Times", 12, L0, nohscroll; fontset = leftheader, "Times", 12, L2; fontset = leftfooter, "Times", 12, L2; fontset = reserved1, "Courier New", 10, L0, nowordwrap, nohscroll;] :[font = title; inactive; preserveAspect; startGroup; nohscroll; center; ] Pseudoprimes and Strong Pseudoprimes :[font = subsection; inactive; preserveAspect; startGroup; Cclosed; nohscroll; ] Pseudoprimes :[font = text; inactive; preserveAspect; nohscroll; ] The tests for primality using pseudoprimes and strong pseudoprimes can be easily implemented in Mathematica. Here are some functions which can be used to find pseudoprimes, strong pseodoprimes and implement a probabilisitic primality test. First, to check if a number is a pseudoprime to some base, we can define :[font = input; preserveAspect; nowordwrap; ] pspQ[a_,n_]:= PowerMod[a,n-1,n] ==1 && !PrimeQ[n] :[font = text; inactive; preserveAspect; nohscroll; ] Note that we are using the built-in primeQ to test if the number is composite. To find the smallest pseudoprime to base 5 which is larger than 5, we can type :[font = input; preserveAspect; startGroup; nowordwrap; ] n=5; While[ !pspQ[5,n], n++]; Print[n] :[font = print; inactive; preserveAspect; endGroup; nowordwrap; ] 124 :[font = text; inactive; preserveAspect; nohscroll; ] We verify if 124 is a pseudoprime to base 5. :[font = input; preserveAspect; startGroup; nowordwrap; ] PowerMod[5,123,124] :[font = output; inactive; formatted; output; preserveAspect; endGroup; nowordwrap; ] 1 ;[o] 1 :[font = text; inactive; preserveAspect; endGroup; nohscroll; ] The code can be modified to find pseudoprimes with different bases and conditions. Regarding Carmichael numbers, it is not practical to verify the definition for a large number. Instead it is better to use the characterization of a Carmichael number in terms of its prime factorization. Naturally, this makes the discovery of large Carmichael numbers more difficult. :[font = subsection; inactive; preserveAspect; startGroup; Cclosed; nohscroll; ] Strong Pseudoprimes :[font = text; inactive; preserveAspect; nohscroll; backColorRed = 65535; backColorGreen = 65535; backColorBlue = 65535; fontColorRed = 0; fontColorGreen = 0; fontColorBlue = 0; plain; fontName = "Times New Roman"; fontSize = 12; ] The strong pseudoprime test for n requires determining the even part of n-1. We want to use this test to find large probable primes, hence we should not incorporate the built-in primality test into the strong pseudoprime test. To eliminate the primes, one should use a program for trial division, which is suitable only up to some limit like 10^{12}. Here is a program which tests if a number n passes the strong pseudoprime test to base a. The program is designed to print the intermediate computations. You may delete the print statements if you like. First, we have to factor n-1= 2^r s. This is the prupose of the first While loop. The program then prints s and r. Then we compute the terms b_0= a^s \pmod{n} and its squares until we reach the last value, a^{n-1}\pmod{n} or we get 1. We first check if a^s\equiv 1\pmod{n}. If this is the case, then no further information can be obtained from the computation and we declare the number to a porbable strong pseudoprime. If this is not the case, the While loop computes succesive squares terminates if one of values is 1 or -1 or we have computed all the values. In the last case, if i=r and we haven't encountered 1 or -1 before then the number must be composite. If i