AOJ 1189 Prime Caves (ICPC国内予選2013問題D)

• 問題文

1から m まで(1 <= m <= 10^6)まで螺旋状に書かれたとき, 数字 n (1 <= n <= m)からスタートして下・右下・左下の移動を繰り返し移動します。そのような経路の中で通った素数の個数が最大の経路を求め, その個数と最後に通った素数を出力します。

ある整数が素数かどうかO(1)で判定できるようにエラトステネスの篩等で最初に素数を列挙しておきます。次に螺旋状に数字を書きます。バグりやすいので気をつけましょう。ここまでできたらDPで解を計算するだけです。
dp[y][x] := (x,y) までたどりつくときに通った素数の個数
is_p[n] := 整数 n が素数のとき 1, そうでないときは 0 を返す
a[y][x] := (x,y) に書かれた数字を返す
として
dp[y+1][x] = max(dp[y][x-1] + is_p[a[y][x-1]], max(dp[y][x] + is_p[a[y][x]], dp[y][x+1] + is_p[a[y][x+1]]))
という漸化式が成り立つと思います。

#include <iostream>
#include <cstdio>
#include <algorithm>
using namespace std;

const int MAX_P = 1001000;
const int MAX_W = 1020;
const int dx[4] = {+1,  0, -1,  0};
const int dy[4] = { 0, -1,  0, +1};

// is_p[n] := 正整数 n が素数かどうか
char is_p[MAX_P];
// a[y][x] := (x,y) の洞穴の番号
int a[MAX_W][MAX_W] = {0};
// dp[y][x] := (x,y) にたどり着くまでの素数洞穴を調べた個数
int dp[MAX_W][MAX_W];
// X[k] := 洞穴の番号 k の x 座標 (Yも同様)
int X[MAX_P];
int Y[MAX_P];

void init(){
    // エラトステネスの篩
    fill(is_p, is_p + MAX_P, 1);
    is_p[0] = is_p[1] = 0;
    for(int i = 2; i < MAX_P; i++){
        if( is_p[i] ){
            for(int j = i + i; j < MAX_P; j += i){
                is_p[j] = 0;
            }
        }
    }
}

void make_douketu(int m){
    for(int y = 0; y < MAX_W; y++){
        for(int x = 0; x < MAX_W; x++){
            a[y][x] = 0;
        }
    }
    int k = 1, x = MAX_W / 2, y = MAX_W / 2, dir = 0, d = 2;
    while( k < MAX_P ){
        for(int i = 0; i < (d / 2); i++){
            a[y][x] = k;
            X[k] = x;
            Y[k] = y;
            x += dx[dir];
            y += dy[dir];
            k++;
            if( m < k ) return;
            if( MAX_W <= x || MAX_W <= y ) break;
        }
        dir = (dir + 1) % 4;
        d++;
        if( MAX_W <= x || MAX_W <= y ) break;
    }
}

int main(){
    init();
    int n, m;
    while( cin >> m >> n, n || m){
        make_douketu(m);
        for(int y = 0; y < MAX_W; y++){
            for(int x = 0; x < MAX_W; x++){
                dp[y][x] = -1;
            }
        }
        int sx = X[n], sy = Y[n];
        dp[sy][sx] = is_p[n];
        
        int ans = 0;
        for(int y = 0; y < MAX_W; y++){
            for(int x = 0; x < MAX_W; x++){
                ans = max(ans, dp[y][x]);
                if( a[y][x] == 0 ) continue;
                
                if( dp[y][x] != -1 ){
                    int my = y + 1;
                    if( MAX_W <= my ) continue;
                    
                    int d[3] = {-1,0,1};
                    for(int i = 0; i < 3; i++){
                        int mx = x + d[i];
                        if( mx < 0 || MAX_W <= mx ) continue;
                        
                        int k = a[my][mx];
                        if( k == 0 ) continue;
                        dp[my][mx] = max(dp[my][mx], dp[y][x] + is_p[k]);
                        ans = max(ans, dp[my][mx]);
                    }
                }
            }
        }
        
        if( ans == 0 ){
            cout << "0 0" << endl;
            continue;
        }
        int ans_n = 0;
        for(int y = 0; y < MAX_W; y++){
            for(int x = 0; x < MAX_W; x++){
                if( a[y][x] == 0 ) continue;
                if( dp[y][x] == ans ){
                    int k = a[y][x];
                    if( is_p[k] ){
                        ans_n = max(ans_n, a[y][x]);
                    }
                }
            }
        }
        cout << ans << " " << ans_n << endl;
    }
}