027. Quadratic primes [소수를 생성하는 2차식의 계수 찾기]

​문제 (English)

Euler discovered the remarkable quadratic formula:

$$ n^{2} + n + 41$$

 

It turns out that the formula will produce 40 primes for the consecutive integer values $0 \leq n \leq 39$. However, when $n=40, 40^{2}+40+41=40(40+1)+41$ is divisible by $41$, and certainly when $n=41, 41^{2}+41+41$ is clearly divisible by 41. The incredible formula $n^{2} - 79n + 1601$ was discovered, which produces 80 primes for the consecutive values $0 \leq n \leq 79$. The product of the coefficients, -79 and 1601, is -126479.

 

Considering quadratics of the form:

$n^{2} +an +b, whaer |a| < 1000 and |b| \leq 1000$

 

where $|n|$ is the modulus/absolute value of $n$

e.g. $|11|=11$ and $|-4|=4$

 

Find the product of the coefficients, $a$ and $b$, for the quadratic expression that produces the maximum number of primes for consecutive values of $n$, starting with $n=0$.

 

Answer: -59231

문제 분석

---

어떤 2차식 $f=n^{2}+an+b$에서 $a$와 $b$를 잘 선택하면, $n=0$에서 f의 값이 소수이고, $n$을 계속 증가해도 어느 순간까지는 계속 소수가 생성된다. 이때 $|a| < 1000$이고 $|b| \leq 1000$인 범위에서, $n$을 계속 증가시킬 때 가장 많은 수의 소수를 연속해서 만들어내는 $a, b$ 값을 찾고, 두 곱 $ab$를 계산하라는 문제이다.

 

$a$와 $b$를 주어진 범위 안에서 증가시키면서 가능한 모든 조합을 찾고, 해당 조합일 때 2차식 $f=n^{2}+an+b$에 대해 $n=0$에서부터 1씩 증가시키면서 $f$가 소수가 안될 때까지 찾으면 되겠다.  

 

풀이 1 

---

$a$와 $b$를 주어진 범위안에서 조합을 만들려면 2중 for 구문을 쓰면 되겠다. 

$a$와  $b$의 범위가 문제에서는 $|a|<1000$이고 $|b| \leq 1000$으로 조금 다르게 되어 있어서, 코드에도 그렇게 했으나, a, b 모두 1000을 포함하게 문제를 내는 게 더 일관성이 있어 보인다. 물론, 그렇게 하더라도 답이 달라지거나 하진 않는다.

 

a, b 조합에 대해 몇 개의 소수가 나오는지는 num_of_primes(a,b) 함수에 의해 구해지고, 여기서 구한 값을 다른 a, b 조합에서 구한 소수 개수와 비교해서, 최댓값을 찾는다.

#[cfg(test)]
fn find_max_coef(k:i32) -> (i32, i32) {
    let mut max = 0;    
    let mut max_a=0; let mut max_b=0;
    for a in -(k-1)..k {
        for b in -(k-1)..=k{
            let t = num_of_primes(a, b);
            if t > max {max=t; max_a=a; max_b=b;}
        }
    }
    return (max_a, max_b);
}

 

num_of_primes(a, b) 함수는 간단하다. 

$n^2+an+b$에서, $n=0$부터 식에 대입해서 나온 값이 소수인지 판가름하면서, 소수가 아닌 값이 나올 때까지 $n$을 증가시킨다. 연속된 소수가 나올 때까지가 원하는 조건이기에, n을 증가시키다가 하나라도 소수가 아닌 값이 나오면, 그 직전까지의 개수가 생성되는 소수의 개수가 된다. 

#[cfg(test)]
fn num_of_primes(a:i32, b:i32) -> usize {
    let mut n:i32 = 0;
    while is_prime(n*n + a*n + b) {n += 1;}
    return n as usize;
}

 

num_of_primes 함수 안에서, 식에 의해 생성된 수가 소수인지 아닌지를 is_prime 함수로 판단한다. 

 

#[cfg(test)]
fn is_prime(n:i32) -> bool {
    if n <= 1 {return false;}    
    if n == 2 {return true;}    
    if n % 2 == 0 {return false;}
    if n % 3 == 0 {return false;}

    let sqrt_n = (n as f64).sqrt() as usize;
    for i in (5..=sqrt_n).step_by(6) {
        if n as usize % i == 0 {return false; }
        if n as usize % (i+2) == 0 {return false;}
    }
    return true;
}

is_prime 함수는 최적화를 꽤 했다. 

2와 3의 배수인지 먼저 검사를 하고, for 루프에서는 $\sqrt{n}$까지만 검사하게 하고, step_by를 6으로 했다. 

for 루프에서 n % i와 n % (i+2)를 하고, for 루프 이전에 2와 3에 대해서 배수인지 확인했기에 step_by(6)이 가능한 것이다.

 

이처럼 각각의 수에 대해서 소수를 판단하지 않고, 에라토스테네스의 체를 이용해서 주어진 수 공간에 대해 소수인지 아닌지를 먼저 판별해 놓고, 어떤 수에 대해 소수여부 판단은 이 체에 의해 한다면 상수 시간에 소수여부를 판단할 수 있다.

그러나, 여기서는 그냥 is_prime 함수를 써서 판단하겠다. 이렇게해도 충분히 속도가 나온다. 

 

풀이 2 

---

풀이1에서 a, b의 조합을 주어진 범위 내 모두에 대해서 수행했다. 

문제에서 주어진 조건을 이용하면, for 루프의 범위를 좀 더 좁힐 수 있다.

 

1) $n=0$에 대해서 $n^2+an+b$가 소수라야하기에, $b$는 소수라야한다. 즉, $b \leq 2$

2) 연속해서 소수라야하기에, 최소한 2개연속 소수이려면 n=1일 때도 소수라야한다. 즉, $1+a+b \leq 2$. 위 1번에서 b의 범위가 좁혀졌기에, b에 의한 a의 범위를 한정해 보면, $a \geq 1-b$

 

위 제약조건에 따라 코드를 짜보면,

#[cfg(test)]
fn find_max_coef1(k:i32) -> (i32, i32) {
    let mut max = 0;    
    let mut max_a=0; let mut max_b=0;

    for b in 2..k {
        for a in (1-b)..k{
            let t = num_of_primes(a, b);
            if t > max {max=t; max_a=a; max_b=b;}
        }
    }
    return (max_a, max_b);
}

 

총평

---

n=0에서 시작해서 n을 증가시키면서 나오는 수에 대해 계속해서 소수가 나와야하기에, 소수가 아닐 때까지만 n을 증가시켜 보면 된다는 것을 인지하는 것이 핵심.

 

코드상으로는, 해당 수가 소수인지 아닌지를 효과적으로 판단하는 것이 핵심

 

전체 소스코드

---

p27.rs 

#[test]
fn test(){    
    // assert_eq!(-41, answer(42));
    assert_eq!(-59231, answer(1000));
}

#[cfg(test)]
//return (a,b) which produce max prime in n^2+an+b, |a|<=k, |b|<=k
fn answer(k:i32) -> i32 {
    let (a,b) = find_max_coef(k);
    return a*b;
}

#[cfg(test)]
fn find_max_coef(k:i32) -> (i32, i32) {
    let mut max = 0;    
    let mut max_a=0; let mut max_b=0;
    for a in -(k-1)..k {
        for b in -(k-1)..=k{
            let t = num_of_primes(a, b);
            if t > max {max=t; max_a=a; max_b=b;}
        }
    }

    println!("max={} a={} b={}",max, max_a, max_b);
    return (max_a, max_b);
}

#[cfg(test)]
fn num_of_primes(a:i32, b:i32) -> usize {
    let mut n:i32 = 0;
    while is_prime(n*n + a*n + b) {n += 1;}
    return n as usize;
}

#[cfg(test)]
fn is_prime(n:i32) -> bool {
    if n <= 1 {return false;}    
    if n == 2 {return true;}    
    if n % 2 == 0 {return false;}
    if n % 3 == 0 {return false;}

    let sqrt_n = (n as f64).sqrt() as usize;
    for i in (5..=sqrt_n).step_by(6) {
        if n as usize % i == 0 {return false; }
        if n as usize % (i+2) == 0 {return false;}
    }
    return true;
}

//--------------------
//최적화
#[test]
fn test1(){    
    assert_eq!(-41, answer1(42));
    assert_eq!(-59231, answer1(1000));
}

#[cfg(test)]
//return (a,b) which produce max prime in n^2+an+b, |a|<=k, |b|<=k
fn answer1(k:i32) -> i32 {
    let (a,b) = find_max_coef1(k);
    return a*b;
}

#[cfg(test)]
fn find_max_coef1(k:i32) -> (i32, i32) {
    let mut max = 0;    
    let mut max_a=0; let mut max_b=0;

    for b in 2..k {
        for a in (1-b)..k{
            let t = num_of_primes(a, b);
            if t > max {max=t; max_a=a; max_b=b;}
        }
    }
    return (max_a, max_b);
}

 

 

hackerrank에 올린 코드는 다음과 같다.

 

use std::io::{self, BufRead};
fn main() {
    let stdin = io::stdin();
    let mut it = stdin.lock().lines();
    let n = it.next().unwrap().unwrap().trim().parse::<i32>().unwrap();
    let (a,b) = answer(n);
    println!("{} {}",a,b);
}

//return (a,b) which produce max prime in n^2+an+b, |a|<=k, |b|<=k
fn answer(k:i32) -> (i32, i32) {
    let (a,b) = find_max_coef(k);
    return (a,b);
}

fn find_max_coef(k:i32) -> (i32, i32) {
    let mut max = 0;    
    let mut max_a=0; let mut max_b=0;
    for a in -(k-1)..k {
        for b in -(k-1)..k{
            let t = num_of_primes(a, b);
            if t > max {max=t; max_a=a; max_b=b;}
        }
    }
    return (max_a, max_b);
}

fn num_of_primes(a:i32, b:i32) -> usize {
    let mut n:i32 = 0;
    while is_prime(n*n + a*n + b) {n += 1;}
    return n as usize;
}

fn is_prime(n:i32) -> bool {
    if n <= 1 {return false;}    
    if n == 2 {return true;}    
    if n % 2 == 0 {return false;}
    if n % 3 == 0 {return false;}

    let sqrt_n = (n as f64).sqrt() as usize;
    for i in (5..=sqrt_n).step_by(6) {
        if n as usize % i == 0 {return false; }
        if n as usize % (i+2) == 0 {return false;}
    }
    return true;
}

끝