Pointwise convergence of bilinear polynomial averages over the primes

We show that on a $\sigma$-finite measure preserving system $X = (X,\nu, T)$, the non-conventional ergodic averages $$\mathbb{E}_{n \in [N]} \Lambda(n) f(T^n x) g(T^{P(n)} x)$$ converge pointwise almost everywhere for $f \in L^{p_1}(X)$, $g \in L^{p_2}(X)$, and $1/p_1 + 1/p_2 \leq 1$, where $P$ is a polynomial with integer coefficients of degree at least $2$. This had previously been established with the von Mangoldt weight $\Lambda$ replaced by the constant weight $1$ by the first and third authors with Mirek, and by the M\"obius weight $\mu$ by the fourth author. The proof is based on combining tools from both of these papers, together with several Gowers norm and polynomial averaging operator estimates on approximants to the von Mangoldt function of ''Cram\'er'' and ''Heath-Brown'' type.