Question

By replacing the 1st digit of *3, it turns out that six of the nine possible values: 13, 23, 43, 53, 73, and 83, are all prime.

By replacing the 3rd and 4th digits of 56**3 with the same digit, this 5-digit number is the first example having seven primes among the ten generated numbers, yielding the family: 56003, 56113, 56333, 56443, 56663, 56773, and 56993. Consequently 56003, being the first member of this family, is the smallest prime with this property.

Find the smallest prime which, by replacing part of the number (not necessarily adjacent digits) with the same digit, is part of an eight prime value family.

Haskell

import Data.List (sort, group)

primes :: [Int]
primes = 2 : sieve primes [3,5..] where
    sieve (p:ps) xs = h ++ sieve ps [x | x <- t, rem x p /= 0]
                      where (h, t) = span (< p*p) xs

isPrime :: Int -> Bool
isPrime n | n < 1 = False
          | otherwise = not $ or [n `rem` x == 0 | x <- [2..floor $ sqrt $ fromIntegral n]]

groupFrequency :: Ord a => [a] -> [(Int, a)]
groupFrequency xs = reverse $ sort [(length cs, head cs) | cs <- group $ sort xs]

candidates :: [(Char, String)]
candidates = [(snd $ head g, s) | p <- primes,
                                  let s = show p,
                                  let g = groupFrequency s,
                                  any ((== 3) . fst) g]

replace :: Eq a => a -> a -> [a] -> [a]
replace old new = foldr (\x acc -> (if x == old then new else x) : acc) []

expand :: (Char, String) -> [String]
expand (d, ds) = [new | n <- ['1'..'9'], let new = replace d n ds, isPrime (read new)]

main :: IO ()
main = putStrLn $ snd . head $ filter ((== 8) . length . expand) candidates

$ ghc -O2 -o prime-replacement prime-replacement.hs
$ time ./prime-replacement
real   0m0.071s
user   0m0.071s
sys    0m0.000s

Python

#!/usr/bin/env python
from itertools import *
def primes(n): 
    if n==2: return [2]
    elif n<2: return []
    s=list(range(3,n+1,2))
    mroot = n ** 0.5
    half=(n+1)//2-1
    i=0
    m=3
    while m <= mroot:
        if s[i]:
            j=(m*m-3)//2
            s[j]=0
            while j<half:
                s[j]=0
                j+=m
        i=i+1
        m=2*i+3
    return [2]+[x for x in s if x]

from bisect import bisect_left
# sqrt(1000000000) = 31622
__primes = primes(31622)
def is_prime(n):
    # if prime is already in the list, just pick it
    if n <= 31622:
        i = bisect_left(__primes, n)
        return i != len(__primes) and __primes[i] == n
    # Divide by each known prime
    limit = int(n ** .5)
    for p in __primes:
        if p > limit: return True
        if n % p == 0: return False
    # fall back on trial division if n > 1 billion
    for f in range(31627, limit, 6): # 31627 is the next prime
        if n % f == 0 or n % (f + 4) == 0:
            return False
    return True

def number_families(num):
    digits = [d for d in str(num)]
    products = list(product((True, False), repeat=len(digits)))[1:-1]
    for p in products:
        pattern = ''
        for i, x in enumerate(p):
            if x:
                pattern += digits[i]
            else:
                pattern += 'X'
        yield [int(pattern.replace('X', str(n))) for n in range(10)]

def main():
    for prime in primes(1000000):
        for number_family in number_families(prime):
            prime_family = [n for n in number_family if is_prime(n) and len(str(n)) == len(str(prime))]
            if len(prime_family) == 8 and prime in prime_family:
                print(prime)
                return

if __name__ == "__main__":
    main()

$ time python3 changing-primes.py
real   0m9.551s
user   0m9.551s
sys    0m0.000s