A pentagonal number is generated by the formula:

\[P_n = n(3n-1)/2\]

Find the pair of pentagonal numbers \(P_j\) and \(P_k\) where the sum \(P_j + P_k\) and \(P_j - P_k\) is also pentagonal. And where \(D = |P_j - P_k|\) is minimised.

The answer we’re looking for is D.

First a function for generating pentagonal numbers:

pent <- function(x){
  x*(3*x - 1)/2
}

Next a nice long list of pentagonal numbers:

pent_numbers <- pent(1:10000)

Third: Generate the pair-wise combinations:

library(tidyr)
pent_comb <- expand_grid(Var1 = pent_numbers, Var2 = pent_numbers)

Fifth: Filter on the criteria:

library(dplyr)
answer <- pent_comb |> 
  filter(abs(Var1 - Var2) %in% pent_numbers) |> 
  filter((Var1 + Var2) %in% pent_numbers) |> 
  filter((Var1 - Var2) %in% pent_numbers) |> 
  summarise(answer = Var1 - Var2) 

A bit lucky that I chose a length of pentagonal numbers that actually included the result in the first attempt.