Are there infinitely many pairs of primes where each divides one more than the square of the other?

I have the following question on number theory that is eating my head.

Are there infinitely many primes $p,q$ such that $p | (q^2 + 1)$ and $q | (p^2 + 1)$?

I can see $13,5$ and $2,5$ has the required property by a simple calculation. A computer code reveals the pair $89,233$ as well.

But I don't know if there are more. Can someone shed some light on this matter?

Thanks!

The part about $p,q$ being prime is a bit of a red herring. I did the same thing without primality required, got this with the larger number up to 46,000. Should not be too difficult to prove that every other Fibonacci number is involved, the ones that are also Markov numbers. http://en.wikipedia.org/wiki/Markov_number From the outcome and context, it is likely that what this site usually calls "Vieta Jumping" is involved, this being fundamental to the Markov numbers as well as being a method for some otherwise intractable contest problems. Essentially, I suggest trying to prove that if you have a solution pair $(p,q)$ with $q$ larger, then there is a solution pair $(q',p)$ with $p$ now the larger. Oh, it is legitimate to call them $p_n,$ and try to prove $p_{n+2} = 3 p_{n+1} - p_n.$ If successful, this will show Fibonacci-ness for all.

     1         2
     2         5
     5        13
    13        34
    34        89
    89       233
   233       610
   610      1597
  1597      4181
  4181     10946
 10946     28657

Alright, my memory was correct.

We have $p | (q^2 + 1)$ and $q | (p^2 + 1),$ with $p,q$ positive integers but not necessarily prime. We get $$ pq | (q^2 + 1)(p^2 + 1) = p^2 q^2 + p^2 + q^2 + 1. $$ Since $pq$ automatically divides $p^2 q^2,$ this means that $$ pq | (p^2 + q^2 + 1). $$ The second example on this Wikipedia page, see See http://en.wikipedia.org/wiki/Vieta_jumping#Example_2 , shows that $$ \color{magenta}{ p^2 - 3 p q + q^2 = -1}. $$ That is actually enough, we can describe all integer solutions of $p^2 - 3pq + q^2 = -1.$ Meanwhile, they are also Markov numbers, because $$ 1 + p^2 + q^2 = 3pq, $$ meaning we have a Markov triple $(1,p,q).$ Indeed, the other two numbers in a Markov triple containing a $1$ are odd-index Fibonacci numbers, see http://en.wikipedia.org/wiki/Markov_number#Markov_tree

EDIT, Saturday morning; There are several ways to prove the claim that the solutions to $p^2 -3pq+p^2 = -1$ involve odd-indexed Fibonacci numbers, as is stated on the Markov page and an OEIS page. I am going to define integer variables $$ x = p - 2q, \; \; y = q $$ so that $$ p^2 - 3pq + q^2 = x^2 + xy - y^2. $$ The diagram below shows the Conway topograph for $x^2 + xy - y^2$ and shows how all solutions for $x^2 + xy - y^2 = -1$ involve Fibonacci numbers, sometimes mutiplied by $(-1).$ The diagram also illustrates the proper automorphism group, which is just two formulas, $$ (x,y) \mapsto (x+y, x+2y) $$ to move to the right, and $$ (x,y) \mapsto (2x-y, -x+y) $$ to move to the left. Anyway, these are all $(x,y)$ pairs that give $(-1),$ and re get back with $$ p = x+2y, \; q = y $$ Conway's little book is available at http://www.maths.ed.ac.uk/~aar/papers/conwaysens.pdf

Oh, Cayley-Hamilton. The automorphism group generator matrix $$ A= \left( \begin{array}{rr} 1 & 1 \\ 1 & 2 \end{array} \right) $$ gives us the rule for the sequence of solution pairs, $$ \left( \begin{array}{rr} 1 & 1 \\ 1 & 2 \end{array} \right) \left( \begin{array}{r} x_n \\ y_n \end{array} \right) = \left( \begin{array}{r} x_{n+1} \\ y_{n+1} \end{array} \right). $$ Do this twice and you get $(x_{n+2}, y_{n+2})^T.$ However, $A$ has determinant $1$ and trace $3,$ so Cayley-Hamilton says $$ A^2 - 3 A + I = 0. $$ In turn, this gives separate rules, $$ x_{n+2} = 3 x_{n+1} - x_n, $$ $$ y_{n+2} = 3 y_{n+1} - y_n. $$ In the opposite direction, actually the same rule $$ x_{n} = 3 x_{n+1} - x_{n+2}, $$ $$ y_{n} = 3 y_{n+1} - y_{n+2}. $$ We also have $$ p = x+2y, \; q = y, $$ therefore we have $$ p_{n+2} = 3 p_{n+1} - p_n, $$ $$ q_{n+2} = 3 q_{n+1} - q_n. $$ $$ p_{n} = 3 p_{n+1} - p_{n+2}, $$ $$ q_{n} = 3 q_{n+1} - q_{n+2}. $$

Finally, Hagen's comments about twin primes, from Wikipedia,

It is not known whether there are infinitely many Fibonacci primes. The first 33 are Fn for the n values (sequence A001605 in OEIS):

3, 4, 5, 7, 11, 13, 17, 23, 29, 43, 47, 83, 131, 137, 359, 431, 433, 449, 509, 569, 571, 2971, 4723, 5387, 9311, 9677, 14431, 25561,

30757, 35999, 37511, 50833, 81839.

In addition to these proven Fibonacci primes, there have been found probable primes for

n = 104911, 130021, 148091, 201107, 397379, 433781, 590041, 593689, 604711, 931517, 1049897, 1285607, 1636007, 1803059, 1968721,

2904353.1

Except for the case n = 4, all Fibonacci primes have a prime index, because if a divides b, then F_a also divides F_b, but not every prime is the index of a Fibonacci prime.

Certainly appears, proof not expected ever, that the prime solutions are the very small ones along with $$ (F_{11}, F_{13}); \; \; (F_{431}, F_{433}); \; \; (F_{569}, F_{571}). $$