32

If not assumed otherwise m and n can be whatever, so you can do e.g. this : Solve[ Prime[n] + Prime[m] == 100, {n, m}, Integers] {{n -> 2, m -> 25}, {n -> 5, m -> 24}, {n -> 7, m -> 23}, {n -> 10, m -> 20}, {n -> 13, m -> 17}, {n -> 15, m -> 16}, {n -> 16, m -> 15}, {n -> 17, m -> 13}, {n -> 20, m -&...


26

Actually, I believe the issue reduced to that of implementing PrimePi[]. It is easy to implement Prime[] using PrimePi[] and FindRoot[] — in fact this is done on page 134 of Bressoud and Wagon, "A Course in Computational Number Theory". So all you need is to have a fast implementation of PrimePi[]. The first efficient way was found by Legendre in 1808. The ...


25

This answer should be read upside down, since the last edit has the fastest, neatest and shortest answer Module[{$guard = True}, happyQ[i_] /; $guard := Block[{$guard = False, appeared}, appeared[_] = False; happyQ[i] ] ] e : happyQ[_?appeared] := e = False; happyQ[1] = True; e : happyQ[i_] := e = (appeared[i] = True; happyQ[#.#&@...


24

PrimeQ and FactorInteger use different algorithms. In general asking whether a number is prime is an easier problem than finding its factors. To quote the documentation, "PrimeQ first tests for divisibility using small primes, then uses the Miller–Rabin strong pseudoprime test base 2 and base 3, and then uses a Lucas test", while "FactorInteger switches ...


23

It's due to an implementation-dependent issue. We should try to improve on it. Has not been much clamor to do so, therefore it has not been a high priority. --- edit --- I've had a look at the code. It is quite intentional that the largest is around what you state (I see the constant being set to $7.783516108362\times 10^{12}$). It has to do with this ...


20

The built-in functionPrimeOmega gives you the number of prime factors and counts multiplicities. Therefore, this can easily be used to give you semi-primes as you have defined them: With[{r = Range[50]}, Pick[r, PrimeOmega[r], 2]]


19

A Trace reveals the problem: Product[If[PrimeQ[p], (p^2 + 1)/(p^2 - 1), 1], {p, 2, Infinity}] // Trace (* {Product[If[PrimeQ[p], (p^2 + 1)/(p^2 - 1), 1], {p, 2, Infinity}], {{PrimeQ[p], False}, If[False, (p^2 + 1)/(p^2 - 1), 1], 1}, Product[1, {p, 2, Infinity}], 1} *) The If[] statement is evaluated before the Product. In turn PrimeQ[p] ...


18

In general Mathematica cannot compute symbolically infinite sums over primes because of the lack of appropriate mathematical tools. However there are infinite products over primes which are basically well understood on the mathematical level. One famous example is the Euler formula for the Riemann zeta function, one of the most beautiful (and mysterious ...


17

I want to answer the part of the question, "How could my son be expected to find a prime factor?" Well, this depends on what your son has been taught, of course. A first thing to notice is that, since 99! is divisible by every prime less than 99, 99! - 1 is not divisible by any of those primes; so 101 is the smallest prime which could be a factor of it. So ...


16

You could use FactorInteger to find out whether or not there are exactly two primes building up a number: SemiPrimeQ[n_Integer] := With[{factors = FactorInteger[n]}, Total[factors[[All, 2]]] == 2 ] The rest is easy: Select[Range[50], SemiPrimeQ] (* {4, 6, 9, 10, 14, 15, 21, 22, 25, 26, 33, 34, 35, 38, 39, 46, 49} *) And for those who like inline ...


16

I believe this is correct, and very fast: fn[x_Integer, n_Integer] := Complement[Range @ x, Join @@ Range[#, x, #]] & @ Prime @ Range @ n Test: fn[10000, 1223] {1, 9929, 9931, 9941, 9949, 9967, 9973} It seems I am a bit late to return to this problem and Simon Woods already provided a memory optimized approach. His sieve is comparatively fast ...


16

This is competitive with Mr Wizards code and seems faster in some cases: fn2[x_Integer, n_Integer] := Module[{y = Range @ x}, (y[[# ;; x ;; #]] = 0) & /@ Prime[Range @ Min[n, PrimePi @ x]]; SparseArray[y]["NonzeroValues"]] AbsoluteTiming[fn[10000, 1223];] (* {0.004000, Null} *) AbsoluteTiming[fn2[10000, 1223];] (* {0.010001, Null} *) ...


14

Turning my comment into an answer, One of the tests performed by PrimeQ for machine-sized integers, namely Miller-Rabin using up to the first 12 primes as bases (as of version 10) has been proved correct for integers up to $2^{64}$ (in fact, the smallest number which that test falsely declares a prime is known to be $3186 65857 83403 11511 67461.$) Of ...


13

This probably counts as cheating, since it uses the fact that all unhappy numbers end up in a cycle including 4. But I like it for simplicity... happyQ[1]=True; happyQ[4]=False; happyQ[n_]:=happyQ[n]=happyQ[#.#&@IntegerDigits[n]] This works with Leonid's happyPrimeN function.


13

If you need all the primes that are twin primes up to n then. TwinPrimesUpTo[n_] := Select[ Most@NestWhileList[ NextPrime, Prime[1], # <= n & ] , (PrimeQ[# + 2]||PrimeQ[# + 2]) &] It was pointed out by @KennyColnago in other answer that it's customary to count the number of pairs of twin primes, in which case a list of the ...


13

out = Select[MatchQ[IntegerDigits[#], {___, 5, ___, 4, ___, 3, ___}] &]@ Prime[Range[50000]]; Length@out 1588 Short@out {1543, 2543, 5413, 5431, 5437, 5443, 5483, 5743, 5843, 8543, <<1568>>, 599243, 599413, 599843, 601543, 602543, 605413, 605443, 605543, 610543, 611543} If you need just the count: count = Count[{___, 5, ___, 4, ___, 3, ...


12

I prefer Bold and Larger in Style: Animate[ Grid[ Partition[ Table[ Style[i, Bold, Larger, If[i > j, Black, If[ PrimeQ @ i, Blue, Gray]]], {i, 100}], 10], Spacings -> {1, 1}], {j, 0, 100}, Paneled -> False] but if you like delete from the ...


12

Simple top-level solution Here is a simplistic completely top-level code: ClearAll[happyQ]; happyQ[n_] := Block[{appeared}, appeared[_] = False; Take[ NestWhileList[ Total[IntegerDigits[#]^2] &, n, (! appeared[#] && (appeared[#] = True)) & ], -2] == {1, 1}]; Clear[happyPrimeN]; ...


11

Clear[happyPrimeN]; happyPrimeN[2000] = 137653; happyPrimeN[2000] // AbsoluteTiming {0., 137653} But seriously, here's a memoised, recursive happyQ that can be used with Leonid's happyPrimeN Clear[sos, happyQ]; sos[k_Integer] := sos[k] = #.# &[IntegerDigits[k]]; happyQ[k_Integer] := happyQ[k] = happyQ[k, {}]; happyQ[1, history_List] := True; happyQ[...


11

I cannot take credit for this, but here is a rather interesting and compact implementation that I found somewhere: F[n_] := {Re[#], Im[#]} & /@ Fold[Join[#1, Last[#1] + I^#2 Range[#2/2]] &, {0}, Range[n]] G[n_] := Table[#[[Prime[k]]], {k, 1, PrimePi[n^2/4 + 1]}] &[F[n]] ListPlot[G[500], AspectRatio -> Automatic, Axes -> None, Frame -> ...


11

In a nutshell, factoring integers is a harder problem than determining primality. This seeming asymmetry is exploited as a component in modern computer security systems (in the form of the RSA cryptosystem). Determining primality has long been known to be doable in polynomial time using a variety of probabilistic algorithms, many of which (as mentioned by ...


11

Take a look at IntegerPartitions, although it relies on brute-force enumeration that is unlikely to scale well. f1[n_] := IntegerPartitions[n, {2}, Prime @ Range @ PrimePi @ n, 1] f2[n_] := Length @ IntegerPartitions[n, {2}, Prime @ Range @ PrimePi @ n] Test: f1[3412] {{3407, 5}} f2[3412] 43 Further experiments Anton Antonov's recent answer ...


11

Prime and IntegerDigits are Listable, ___ (i.e. BlankNullSequence) is a pattern object that can stand for any sequence of zero or more expressions. This is a nice approach and quite fast: Length @ Cases[ IntegerDigits @ Prime @ Range @ 50000, {___, 5, ___, 4, ___, 3, ___}] 1588 Remarks This metod and kglr's operator approach ...


10

Here's my take: (* Brent's algorithm for cycle detection *) happyQ[start_Integer] := Module[{cyc, f, hare, pow, tortoise}, f = Total[IntegerDigits[#]^2] &; cyc = pow = 1; tortoise = start; hare = f[start]; While[tortoise =!= hare, If[pow == cyc, tortoise = hare; pow *= 2; cyc = 0;]; ...


10

Might be faster to do batches by checking for nontrivial gcds. Below I group 100 primes into each batch, take their products, and do a gcd computation until finding that there is a divisor. The example sieves a million primes. random = RandomPrime[10^1000]; list = Table[k random + 1, {k, 1, 2000}]; pprods = Apply[Times, Partition[Prime[Range[10^6]], 10^3], ...


10

Artes's answer is just fine. This variation works as follows. When you click on a number, the background of that number turns yellow and that of each of its divisors turns light blue. DynamicModule[{s = 101}, Grid[Partition[Dynamic@Button[Style[#, 16], (s = #), Background -> Which[ # == s, Yellow, Divisible[s, #], LightBlue, True, White], Appearance -...


10

A naive approach would be this: primePower[n_] := Count[ Range @ n, _?PrimePowerQ] This function works well however it might be very inefficient for large n. It takes a bit to evaluate e.g. primePower[10^6] 78734 which is only a little bigger than PrimePi[10^6] 78498 The latter is much more efficient since it uses advanced algorithms for counting ...


10

There are several problems with your code. The first one is that you are missing a couple of semicolons to suppress output and delineate substatements in a compound function. The second problem is that you are trying to assign a new value to x within the function definition. This doesn't work. x already has the value of whatever number you give. You need ...


10

The problem with your code is that N[Product[...]] calls NProduct[...] which turns integers into reals before Prime can evaluate, causing the error: See Details and Options section of documentation of NProduct. For your calculation, you actually do not need these. It looks like the result converges around $0.660162$: Table[N[Product[(Prime[i] (Prime[i] - 2)...


10

You are right!! Reverse /@ Table[IntegerDigits[Prime[i]], {i, 549, 2, -1}] // Flatten // FromDigits // PrimeQ True 7693749334931393929332939193719311937093988318837783368335831583748333\ 8332831283308379733973977396737673167393733373727391739073107379631963\ 7763376317639563346373631363326371633163706339533853185317539553755374\ ...


Only top voted, non community-wiki answers of a minimum length are eligible