A curious prime fraction formula

Question

Given a positive integer \$n\$ output the integers \$a\$ and \$b\$ (forming reduced fraction \$a/b\$) such that:

$$\frac a b = \prod ^n _{k=1} \frac {p^2_k - 1} {p^2_k + 1}$$

Where \$p_k\$ is the \$k\$ th prime number (with \$p_1 = 2\$).

Examples:

1   -> 3, 5
2   -> 12, 25
3   -> 144, 325
4   -> 3456, 8125
5   -> 41472, 99125
15  -> 4506715396450638759507001344, 11179755611058498955501765625
420 -> very long

Probabilistic prime checks are allowed, and it's ok if your answer fails due to limitations in your language's integer type.

---

Shortest code in bytes wins.

Answer 1

M, 9 bytes

RÆN²‘İḤCP

Try it online!

Trivia

Meet M!

M is a fork of Jelly, aimed at mathematical challenges. The core difference between Jelly and M is that M uses infinite precision for all internal calculations, representing results symbolically. Once M is more mature, Jelly will gradually become more multi-purpose and less math-oriented.

M is very much work in progress (full of bugs, and not really that different from Jelly right now), but it works like a charm for this challenge and I just couldn't resist.

How it works

RÆN²‘İḤCP  Main link. Argument: n

R          Range; yield [1, ..., n].
 ÆN        Compute the kth primes for each k in that range.
   ²‘      Square and increment each prime p.
     İ     Invert; turn p² + 1 into the fraction 1 / (p² + 1).
      Ḥ    Double; yield 2 / (p² + 1).
       C   Complement; yield 1 - 2 / (p² + 1).
        P  Product; multiply all generated differences.

Answer 2

Mathematica, 32 bytes

1##&@@(1-2/(Prime@Range@#^2+1))&

An unnamed function that takes integer input and returns the actual fraction.

This uses the fact that (p2-1)/(p2+1) = 1-2/(p2+1). The code is then golfed thanks to the fact that Mathematica threads all basic arithmetic over lists. So we first create a list {1, 2, ..., n}, then retrieve all those primes and plug that list into the above expression. This gives us a list of all the factors. Finally, we multiply everything together by applying Times to the list, which can be golfed to 1##&.

Alternatively, we can use Array for the same byte count:

1##&@@(1-2/(Prime~Array~#^2+1))&

Answer 3

Python 2, 106 bytes

from fractions import*
n=input()
F=k=P=1
while n:b=P%k>0;n-=b;F*=1-Fraction(2*b,k*k+1);P*=k*k;k+=1
print F

The first and fourth lines hurt so much... it just turned out that using Fraction was better than multiplying separately and using gcd, even in Python 3.5+ where gcd resides in math.

Prime generation adapted from @xnor's answer here, which uses Wilson's theorem.

Answer 4

Ruby, 122 77 65 bytes

Thanks to Sherlock for shaving off 10 bytes.

require'prime'
->n{Prime.take(n).map{|x|1-2r/(x*x+1)}.reduce(:*)}

Defines an anonymous function that takes a number and returns a Rational.

Answer 5

PARI/GP, 33 bytes

n->prod(i=1,n,1-2/(prime(i)^2+1))

Alternate version (46 bytes):

n->t=1;forprime(p=2,prime(n),t*=1-2/(p^2+1));t

Non-competing version giving the floating-point (t_REAL) result (38 bytes):

n->prodeuler(p=2,prime(n),1-2/(p^2+1))

Answer 6

Jelly, 14 13 bytes

RÆN²µ’ż‘Pµ÷g/

Try it online! Thanks to @Dennis for -1 byte.

R                       Range [1..n]
 ÆN                     Nth prime
   ²                    Square
    µ                   Start new monadic chain
     ’ż‘                Turn each p^2 into [p^2-1, p^2+1]
        P               Product
         µ              Start new monadic chain
          ÷             Divide by...
           g/           Reduce GCD

Answer 7

Pyth, 26 25

/RiFN=N*MCm,tdhd^R2.fP_ZQ

Try it here or run the Test Suite.

1 byte saved thanks to Jakube!

Pretty naive implementation of the specifications. Uses the spiffy "new" (I have no idea when this was added, but I've never seen it before) P<neg> which returns whether the positive value of a negative number is prime or not. Some of the mapping, etc can probably be golfed...

Answer 8

Julia, 59 42 bytes

n->prod(1-big(2).//-~primes(2n^2)[1:n].^2)

This is an anonymous function that accepts an integer and returns a Rational with BigInt numerator and denominator.

We begin by generating the list of prime numbers less than 2n² and selecting the first n elements. This works because the nth prime is always less than n² for all n > 1. (See here.)

For each p of the n primes selected, we square p using elementwise power (.^2), and construct the rational 2 / (p + 1), where 2 is first converted to a BigInt to ensure sufficient precision. We subtract this from 1, take the product of the resulting array of rationals, and return the resulting rational.

Example usage:

julia> f = n->prod(1-big(2).//-~primes(2n^2)[1:n].^2)
(anonymous function)

julia> f(15)
4506715396450638759507001344//11179755611058498955501765625

Saved 17 thanks to Sp3000!

Answer 9

Convex, 28 bytes

Convex is a new language that I am developing that is heavily based on CJam and Golfscript. The interpreter and IDE can be found here. Input is an integer into the command line arguments. Indexes are one-based. Uses the CP-1252 encoding.

,:)_{µ²1-}%×\{µ²1+}%׶_:Ðf/p

You may or may not consider this answer to be competing since I was working on a few features that this program uses before the challenge was posted, but the commit was made once I saw this challenge go out.

Answer 10

MATL, 18 bytes

:Yq2^tqpwQpZd1Mhw/

Try it online!

Fails for large inputs because only integers up to 2^52 can be accurately represented internally.

Explanation

:     % implicitly take input n. Generate range [1,...,n]
Yq    % first n prime numbers
2^    % square
tqp   % duplicate. Subtract 1. Product
wQp   % swap. Add 1. Product
Zd    % gcd of both products
1M    % push the two products again
h     % concatenate horizontally
w/    % swap. Divide by previously computed gcd. Implicitly display

Answer 11

Mathematica, 45 bytes

Times@@Array[(Prime@#^2-1)/(Prime@#^2+1)&,#]&

Primes? Fractions? Mathematica.

Answer 12

Haskell, 53 bytes

Anonymous function, 53 characters:

(scanl(*)1[1-2%(p*p+1)|p<-nubBy(((>1).).gcd)[2..]]!!)

Try it here (note: in standard GHCi you need first to make sure Data.Ratio and Data.List are imported):

λ (scanl(*)1[1-2%(p*p+1)|p<-nubBy(((>1).).gcd)[2..]]!!) 5
41472 % 99125
:: Integral a => Ratio a

Haskell's list indexing !! is 0-based. (___!!) is an operator section, forming an anonymous function so that (xs !!) n == xs !! n.

It's four bytes less to generate the whole sequence:

λ mapM_ print $ take 10 $     -- just for a nicer output
    scanl(*)1[1-2%(n*n+1)|n<-[2..],all((>0).rem n)[2..n-1]]
1 % 1
3 % 5
12 % 25
144 % 325
3456 % 8125
41472 % 99125
3483648 % 8425625
501645312 % 1221715625
18059231232 % 44226105625
4767637045248 % 11719917990625
:: IO ()

Answer 13

Japt v2.0a0, 25 bytes

_j}jU £[X²ÉX²Ä]ÃyÈ×ÃË÷Fry

Try it

Answer 14

Vyxal, 11 bytes

ʁǎ²₍‹›vΠ:ġḭ  # main program

ʁ            # range 0 to input
 ǎ²          # ith prime, sqaured
   ₍‹›       # push the increment and decrement of that list wrapped
      vΠ     # take the product of both lists
        :ġ   # duplicate, get the gcd
          ḭ  # divide both numbers by the gcd, print implicitly

I don't think the strategy of 1-2/p^2-1 works in Vyxal because it would just end with a floating point number rather than a fraction, but feel free to prove me wrong. (If it did work, it would be a literal translation from Dennis' M answer, but with a byte shaved off for the nth primes, coming out to 8)

Try it Online!

Answer 15

Seriously, 25 bytes

,r`PªD;⌐k`M┬`π`Mi│g;)@\)\

Outputs a\nb (\n is a newline). Large inputs will take a long time (and might fail due to running out of memory) because prime generation is pretty slow.

Try it online!

Explanation:

,r`PªD;⌐k`M┬`π`Mi│g;)@\)\
,r                         push range(input)
  `PªD;⌐k`M                map:
   P                         k'th prime
    ª                        square
     D                       decrement
      ;                      dupe
       ⌐                     add 2 (results in P_k + 1)
        k                    push to list
           ┬               transpose
            `π`M           map product
                i│         flatten, duplicate stack
                  g;)      push two copies of gcd, move one to bottom of stack
                     @\    reduce denominator
                       )\  reduce numerator

Answer 16

Husk, 11 bytes

Πmo§/→←□↑İp

Try it online!

Answer 17

J, 18 bytes

[:*/1-2%1+*:@p:@i.

Try it online!

Uses the same formula that Martin Ender's answer does.