# Tag Info

19

There are many ways to proceed, the best one uses FrobeniusSolve : I Since we know, that a x + b == y /. Solve[{-4 a + b == 11, 16 a + b == -1}, {a, b}] // Simplify {3 x + 5 y == 43} we find FrobeniusSolve[ {3, 5}, 43] {{1, 8}, {6, 5}, {11, 2}} a bit more straightforward way : II {x, y} /. Solve[ (a x + b == y /. Solve[ {-4 a + b == 11, 16 ...

13

This is caused by a bug in RootReduce for Root objects representing last coordinates of solutions of triangular systems. The bug affects cases where the last coordinate of the solution is real, but some of the other coordinates are not real. Thanks for pointing it out. The problem can be fixed with the following patch (you can put it in your init.m file). ...

13

Is this what you are searching for? a = {-4, 11}; b = {16, -1}; dy = (b[[2]] - a[[2]])/(b[[1]] - a[[1]]); offset = u /. Solve[a[[2]] == dy*a[[1]] + u, u][[1]]; coords = {x, y} /. {Reduce[y == dy*x + offset && x > 0 && y > 0, {x, y}, Integers] // ToRules} (* {{1, 8}, {6, 5}, {11, 2}} *) Graphics[{PointSize[Large], ...

13

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 ...

11

You can also use InterpolatingPolynomial with Solve, Reduce or Eliminate: a = {-4, 11}; b = {16, -1}; coords = Solve[y == InterpolatingPolynomial[{a, b}, x] && 0 <= x <= 16&&0<=y, {x, y}, Integers][[All, All, 2]]; (* or *) coords={ToRules[Reduce[ y == InterpolatingPolynomial[{a, b}, x] && 0 <= x <= ...

10

It's meant to be done divide-and-conquer style. Here is one way to go about that. listMod[n_, {val_}] := {Mod[n, val]} listMod[n_, vals : {_, __}] := With[{len = Length[vals], rem = Mod[n, Times @@ vals]}, Join[listMod[rem, Take[vals, Floor[len/2]]], listMod[rem, Drop[vals, Floor[len/2]]]] ] Your example: n = 31415926535; primeslist = ...

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 ...

8

(nextPrime[#1] = #2) & @@@ {{-3, 2}, {-2, 2}, {-1, 2}, {0, 2}, {1, 2}, {2, 3}}; nextPrime[n_Integer?EvenQ] := nextPrime[n - 1]; nextPrime[n_Integer] /; PrimeQ[n + 2] := n + 2; nextPrime[n_Integer] := nextPrime[n + 2] nextPrime[n_ /; n \[Element] Reals] := nextPrime[Floor@n]

8

Even though Mathematica has a broad range of powerful capabilites (see e.g. this comparison of computer algebra systems) in related fields (number theory, quantifier elimination) it sometimes doesn't appear to be clever enough to prove simple theorems, e.g. this should yield False however we get back the input: Resolve[ Exists[p, p ∈ Integers && ...

8

Number theory questions are always a huge accumulator for up votes. :) From my experience I can say that the builtin MangoldtLambda function is pretty slow. So let's define a Mangoldt function on our own. The Mangoldt function is defined by: $\Lambda(n) \equiv \left\{ \begin{array}{1 1} ln\ p & \quad \text{if n =$p^k$for p a prime}\\ 0 & \quad ... 7 Following Daniel's comment, I post my solution: PolynomialCRT[pol_List, mod_List, x_] := Module[{m, q, l}, l = Length[pol]; m = Table[Times @@ Drop[mod, {i}], {i, l}]; q = Table[PolynomialExtendedGCD[m[[i]], mod[[i]], x][[2, 1]], {i, l}]; Simplify[(q pol).m] ] There is no error checking code (for instance, pol and mod have the same length). 7 Perfect numbers: Select[ Range[10^6], Total[Divisors @ #] == 2 # &] {6, 28, 496, 8128} abundant numbers: Select[ Range[10^3], Total[ Most @ Divisors @ #] > # &]//Short { 12, 18, 20, 24, 30, 36, 40, 42, 48, <<228>>, 968, 972, 978, 980, 984, 990, 992, 996, 1000} I used Short to to get only a few since there are: ... 6 I had a clever idea for how to do this using LatticeReduce[], but I decided to code up the Smith-Cornacchia algorithm first to bench mark against, and it was effectively instant for inputs in your range. Here is a sloppy implementation. In particular, I am embarrassed by applying Divisors[] to something which is computed as a product. However, the result is ... 6 I wrote this answer as I was figuring things out. If you just want the answer, copy the definitions of fundQ, Ast, ksFactors, foo, makeOneIndex and magicJ out of the code blocks below. You should have DirichletCharacter[d, magicJ[d], n] == KroneckerSymbol[d, n] Disclaimer 1: This answer is based on plausible interpretations of things not quite said in the ... 6 I get this result as well in both V8 and V9. Product[n^MoebiusMu[n], {n, 1, Infinity}] (* Out: 1/(4*Pi^2) *) It's a simple fact, though, that an infinite product can converge to a non-zero value only if the general term tends to 1. As MoebiusMu takes each of the values$\pm 1(as well as zero) infinitly often, this product simply can't converge. We ... 6 I've found that NSum[] takes a bit too long here to compute Riemann's prime-counting function, so I've resorted to generating the terms and summing them: With[{n = 8*^3}, Total[MoebiusMu[Range[n]] N[LogIntegral[1000^(1/Range[n])]/Range[n]], Method -> "CompensatedSummation"]] 168.35915686601484 which gives a result close to that of ... 6 We can use Reduce to give us: expr = (2^p - (2^2) (3^2))/(3^3); Reduce[{expr == n, n ∈ Integers}, p, Integers] (n | p) ∈ Integers && n >= -1 && p == Log2[9 (4 + 3 n)] Since the only x such that Log2[x] ∈ Integers are powers of two, we must have 9 (4 + 3 n) equal a power of two while n is simultaneously an integer. This clearly ... 6 A faster approach to finding Perfect numbers using DivisorSigma Select[Range[10^6], DivisorSigma[1, #] == 2 # &] {6, 28, 496, 8128} Here's an even faster approach: Pick[#, MapThread[Equal, {DivisorSigma[1, #], 2 #}], True] &[Range[10^6]] and a little bit faster: Pick[#, DivisorSigma[1, #] - 2 #, 0] &@Range[10^6] For Abundant ... 6 Here's a fancy memoized solution: Clear[primePowerCount, primePowerCountcache, iPrimePowerCount] iPrimePowerCount[n1_, n2_] := Count[Range[n1, n2], _?PrimePowerQ] primePowerCountcache = {1}; primePowerCount[1] = 0; primePowerCount[n_?Positive] := Module[{n0, res}, n0 = First@Nearest[primePowerCountcache, n]; If[n0 < n, res = ... 5 Not too hard; all that's needed is a simple application of matrix identities: ColumnHermiteDecomposition[mat_ /; MatrixQ[mat, IntegerQ]] := Transpose /@ HermiteDecomposition[Transpose[mat]] Test: mat = {{1, 2, 3, 2, 2}, {1, 2, 3, 4, 0}, {0, 5, 4, 2, 1}, {3, 2, 4, 0, 2}}; {u, t} = ColumnHermiteDecomposition[mat]; u {{8, 24, 22, ... 5 Without considering the algorithm the only thing that stands out to me is the use of Delete/Partition/Range rather than the native Drop function, as specific native functions are usually faster. A complication of using Drop is that you will need a different exit condition for the loop. That would look something like this: lucky[z_List] := Module[{i = ... 5 (This is a math question, not a Mathematica question.) To add to Artes's answer, there is the well-known(!) identity $$\zeta(-n)=\frac{(-1)^n}{n+1}B_{n+1}$$ so you might as well ask why \begin{align*} -\frac12\times B_2&=-\frac1{14}\times B_{14}\\ -\frac12\times \frac16&=-\frac1{14}\times \frac76 \end{align*} A justification for the ... 5 In order to understand the issue, we should provide the underlying definitions. Mathematica helps in verifying appropriate relations and definitions. The main functional equation relating Riemann's zeta function\zeta\;$, to Euler's$\Gamma\;$, established in Riemann's famous paper Über die Anzahl der Primzahlen unter einer gegebener Grösse (1859, English ... 5 In short NSum cannot handle this sort of sequence. Indeed, strictly the series is not convergent, and some notion of summability/regularization needs to be chosen. Given the nature of MoebiusMu, "Dirichlet" seems appropriate: Sum[MoebiusMu[k], {k, 1, \[Infinity]}, Regularization -> "Dirichlet"] (* -2 *) Here's how one can see NSum is not working ... 5 Here's a functional way to use the property (the property, which has been removed from the original question, was$N'/D' = N/D + 1/D'D$or equivalently$N'D-D'N=1$): farey1[n_] := NestWhileList[ With[{num0 = Numerator[#], den0 = Denominator[#]}, First @ Minimize[{num/den, num den0 - num0 den == 1 && 1 <= den <= n && 1 ... 5 Here's a way to exploit the mediant property of the Farey series. To calculate the mediant: med[{a_, b_}] := (Numerator[a] + Numerator[b])/(Denominator[a] + Denominator[b]); Then the Farey series is: farey[n_] := farey[n] = DeleteCases[ Riffle[ farey[n - 1], med /@ Partition[farey[n - 1], 2, 1]], _?(Denominator[#] > n &)]; with initial ... 5 If I understand correctly: f[n_] := 13 + Ceiling[n, 4]/2; f[Range[20]] {15, 15, 15, 15, 17, 17, 17, 17, 19, 19, 19, 19, 21, 21, 21, 21, 23, 23, 23, 23} More general approach: sample = {15, 15, 15, 15, 17, 17, 17, 17, 19, 19, 19, 19, 21, 21, 21, 21}; linrec = FindLinearRecurrence[sample] {1, 0, 0, 1, -1} f2[n_] := LinearRecurrence[linrec, ... 5 I can't complete with Artes's mathematical knowledge and approach, but simply as a point of reference, for formulating a brute-force approach it will be more memory efficient to use Sum, though it will still be very slow for large input. Sum[Boole @ PrimePowerQ @ i, {i, 5*^6}] // Timing MaxMemoryUsed[] {46.535, 348940} 15083688 4 Here's Eric Weisstein's implementation from MathWorld: Primorial[0] := 1; Primorial[1] := 2; Primorial[n_] := Primorial[n] = Prime[n] Primorial[n - 1]; ChebyshevTheta[n_] := Log[Primorial[PrimePi[n]]] Plot[ChebyshevTheta[x], {x, 0, 100}] Many MathWorld pages have attached notebooks linked near the top of the page. That's where this came from. 4 For reference, here is the v7 code behind NextPrime, which is hard to read before stripping all the private context names. NextPrime[1]; (* preload the definition *) Unprotect[NextPrime]; ClearAttributes[NextPrime, ReadProtected];$Context = "NumberTheoryNextPrimeDump`"; FullDefinition[NextPrime] Yields: Attributes[NextPrime] = {Listable} NextPrime[-3] ...

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