2019年11月2日土曜日

最長共通部分列問題とその解法アルゴリズムとc++のコード

最長共通部分列問題の例題とその解法アルゴリズムについてまとめます。

最長共通部分列問題の例題

AIZU ONLINE JUDGEから引用します

最長共通部分列

最長共通部分列問題 (Longest Common Subsequence problem: LCS)は、
2つの与えられた列 X={x1,x2,...,xm}と Y={y1,y2,...,yn}の最長共通部分列を求める問題です。

ある列 Zが X と Y 両方の部分列であるとき、Z を X とY の共通部分列と言います。
例えば、X={a,b,c,b,d,a,b}, Y={b,d,c,a,b,a} とすると、列 {b,c,a} は X と Y の共通部分列です。
一方、列 {b,c,a} は X と Y の最長共通部分列ではありません。

なぜなら、その長さは 3 であり、長さ 4 の共通部分列 {b,c,b,a} が存在するからです。
長さが 5 以上の共通部分列が存在しないので、列 {b,c,b,a} は X と Yの最長共通部分列の1つです。

与えられた2つの文字列 X、Yに対して、最長共通部分列 Zの長さを出力するプログラムを作成してください。
与えられる文字列は英文字のみで構成されています。

入力

複数のデータセットが与えられます。最初の行にデータセットの数 q が与えられます。続く 2×q 行にデータセットが与えられます。
各データセットでは2つの文字列 X, Yがそれぞれ1行に与えられます。

出力

各データセットについて X, Y の最長共通部分列 Zの長さを1行に出力してください。

制約

1 ≤ q ≤ 150

1 ≤ X,Yの長さ ≤ 1,000

Xまたは Y の長さが 100 を超えるデータセットが含まれる場合、qは 20 以下である。

入力例1

3
abcbdab
bdcaba
abc
abc
abc
bc

出力例1

4
3
2

最長共通部分列(LSC)問題の考え方

2つの文字列$\{x_1,x_2...,x_i\}$を$X_i$、$\{y_1,y_2...,y_j\}$を$Y_j$と表すことにします。
サイズがそれぞれm,nの2つの列X,YのLSCは$X_m$と$Y_n$のLCSを求めることで得られます。
これを部分問題に分割して考えます。

$X_mとY_n$のLCSを求める時は、以下の2つの場合を考えます。

$X_m = y_n$の場合

$X_m = y_n$の場合は、$X_{m - 1}$と$Y_{n - 1}$のLCSに$x_m(= y_n)$を連結したものが$X_m$と
$Y_n$のLCSとなります。

例えば、$X = \{a,b,c,c,d,a\},Y = \{a,b,c,b,a\}$のとき$x_m = y_n$なので、
$X_{m - 1}$と$Y_{n - 1}$のLCSである{a,b,c}に$X_m(= a)$を連結したものが$X_m$と$Y_n$のLCSとなる。

$x_m \neq y_n$の場合

$x_m \neq y_n$の場合は、$X_{m - 1}$と$Y_n$のLCSあるいは$X_m$と$Y_{n - 1}$のLCSのどちらか長いほうが、 $X_m$と$Y_n$のLCSになります。

例えば、$X = \{a,b,c,c,d,b\},Y = \{a,b,c,b,a\}$のとき$X_{m - 1}$と$Y_n$のLCSは${a,b,c}$、
$X_m$と$Y_{n - 1}$のLCSは$\{a,b,c,b\}$なので、$X_m$と$Y_{n - 1}$のLCSが$X_m$と$Y_n$のLCSになります。

最長共通部分列(LSC)解法アルゴリズム

次のような変数を用意して、LCSの部分問題の解を求めていきます

c[m + 1][n + 1] c[i][j]を$X_i$と$Y_j$のLCSの長さとする2次元配列

$c[i][j]$の値は以下のケースの漸化式で求めることができます。

$c[i][j] = 0$ {if i == 0 or j == 0}

$c[i][j] = c[i - 1][j - 1] + 1 $ {if i,j > 0 and $x_i = y_j$}

$c[i][j] = max(c[i][j - 1],c[i - 1][j])$ {if i,j > 0 and $x_i \neq y_j$}

最長共通部分列(LSC)解法疑似アルゴリズム

これらに基づいて、2つの列XとYのLCSを動的計画法で求める疑似アルゴリズムは以下のようになります

lsc(X,Y)
 m = X.length
 n = Y.length

 for i = 1 to m
   c[i][0] = 0
 for j = 1 to n
   c[0][j] = 0
 for i = 1 to m
   for j = 1 to n
     if X[i] == Y[j]
       c[i][j] = c[i - 1][j - 1] + 1
     else if c[i - 1][j] >= c[i][j - 1]
       c[i][j] = c[i - 1][j]
     else
       c[i][j] = c[i][j - 1]

最長共通部分列(LSC)アルゴリズムの計算量

nとmの2重ループより、$O(mn)$のアルゴリズムになります。

最長共通部分列問題のc++による回答

最後に例題の回答を書きます。


#include<iostream>
#include<algorithm>
#include<string>

using namespace std;

#define N 1000

int lcs(string x,string y){
    int c[N + 1][N + 1] = {};
    int maxl = 0;
    x = ' ' + x; // X[0]に空白を挿入
    y = ' ' + y; // Y[0]に空白を挿入

    for(int i = 1;i <= x.size();++i){
        for(int j = 1;j <= y.size();++j){
            if(x[i] == y[j]){
                c[i][j] = c[i - 1][j - 1] + 1;
            }else{
                c[i][j] = max(c[i - 1][j],c[i][j - 1]);
            }
            maxl = max(maxl,c[i][j]);
        }
    }
    
    for(int i = 1;i <= x.size();++i){
        for(int j = 1;j <= y.size();++j){
            if(x[i] == y[j]){
                c[i][j] = c[i - 1][j - 1] + 1;
            }else{
                c[i][j] = max(c[i - 1][j],c[i][j - 1]);
            }
            maxl = max(maxl,c[i][j]);
        }
    }
    // 空白一致分
    return maxl - 1;
}

int main(){
    string s1,s2;
    int n;
    cin >> n;
    
    for(int i = 0;i < n;++i){
        cin >> s1 >> s2;
        cout << lcs(s1,s2) << endl;
    }
}

2019年10月29日火曜日

連鎖行列積問題の考え方

連鎖行列積問題の解法の考え方についてまとめます

例題として、AIZU ONLINE JUDGEの問題を引用します。

問題

n個の行列の連鎖 M1,M2,M3,...,Mn が与えられたとき、
スカラー乗算の回数が最小になるように積 M1M2M3...Mnの計算順序を決定する問題
を連鎖行列積問題(Matrix-Chain Multiplication problem)と言います。

n個の行列について、行列 Mi の次元が与えられたとき、積 M1M2...Mn
の計算に必要なスカラー乗算の最小の回数を求めるプログラムを作成してください。

入力

入力の最初の行に、行列の数 nが与えられます。
続く n 行で行列 Mi(i=1...n) の次元が空白区切りの2つの整数 r、c で与えられます。
r は行列の行の数、cは行列の列の数を表します。

出力

最小の回数を1行に出力してください。

行列の積の考え方

下図のようなl × mの行列Aとm × nの行列Bの積はl × nとなり、
各要素の$C_{ij}$は次の式で得られる

$C_{ij} = \displaystyle \sum_{ i = 1 }^{ m } a_{ik}b_{kj}$

なので、この計算では、l × m × n回のかけ算が必要になる。

複数行列積について考える

下図のような$(M_1,M_2,M_3..M_6)$のような、$M_i$が$p_{i-1}行$$p_i$行の行列であるような、
n個の行列のかけ算を考えます。

これらの行列は様々な順番で計算できるので、行列同士のかけ算の順番によって、
かけ算の計算回数が変わってきます。

連鎖行列問題の解き方考察

連鎖行列において、全ての可能な順番をを調べる方法は$0(n!)$になります。
ただ、小さな部分問題として分割できるので、動的計画法を使うことができます。

動的計画法で解くための考え方

$(M_1,M_2)$を計算するための方法は1通りであり、$p_0 × p_1 × p_2$回のかけ算が必要です。
同様に$(M_2,M_3)$を計算する方法も1通りで、$p_1 × p_2 × p_3$回のかけ算が必要です。
これらのかけ算をコストとして、記録しておきます。

次に、$(M_1,M_2,M_3)$、$(M_2,M_3,M_4)$、...$(M_{n - 2},M_{n - 1},M_n)$ を計算するための最小計算回数つまり最適解を求めます。

例として、$(M_1,M_2,M_3)$の最適解を求めるには、$(M_1(M_2×M_3))$と$((M_1×M_2)M_3))$のコストを計算し、
そのうちの小さい方のコストを$(M_1,M_2,M_3)$のコストとして記録します。

計算量を考えると、

$(M_1(M_2×M_3))$のコスト = $(M_1)$のコスト + $(M_2M_3)$のコスト + $p_0 × p_1 × p_3$
$((M_1×M_2)M_3))$のコスト = $(M_1M_2)$のコスト + $(M_3)$のコスト + $p_0 × p_2 × p_3$

となります。

上記の計算に関して、$(M_1M_2)$のコストと$(M_2M_3)$のコストは記録してあるので、再計算する必要がありません。

また、$1 \leqq i \leqq n$について、$(M_i)$のコストは0になります。

連鎖行列の最適解

一般に連鎖行列の最適解$(M_iM_{i + 1}...M_j)$の最適解は、
$i \leqq k < j$における$(M_iM_{i + 1}...M_k)(M_{k + 1}...M_j)$の最小コストになります

例をあげると
$(M_1M_2M_3M_4M_5)(i = 1,j = 5)$の最適解は

$(M_1)(M_2M_3M_4M_5)$のコスト = $(M_1)$のコスト + $(M_2M_3M_4M_5)$のコスト + $p_0 × p_1 × p_5$(k = 1)のとき

$(M_1M_2)(M_3M_4M_5)$のコスト = $(M_1M_2)$のコスト + $(M_3M_4M_5)$のコスト + $p_0 × p_2 × p_5$(k = 2)のとき

$(M_1M_2M_3)(M_4M_5)$のコスト = $(M_1M_2M_3)$のコスト + $(M_4M_5)$のコスト + $p_0 × p_3 × p_5$(k = 3)のとき

$(M_1M_2M_3M_4)(M_5)$のコスト = $(M_1M_2M_3M_4)$のコスト + $(M_5)$のコスト + $p_0 × p_4 × p_5$(k = 4)のとき

のうちの最小値になります。

上記のアルゴリズムの実装例

以下のような役割の変数を定義します。

m[n + 1][n + 1] m[i][j]を$(M_iM_{i + 1}...M_j)$を計算するための
最小のかけ算の回数とする2次元配列
p[n + 1] $M_i$がp[i - 1] × p[i]の行列となるような1次元配列

これらの変数を用いて、m[i][j]の値を以下のように設定します。

if(i == j){
  m[i][j] = 0 
}

if(i < j){
 m[i][j] = $min_{i \leqq k < j}\{m[i][k] + m[k + 1][j] + p[i - 1] × p[k] × p[j]\}$
}

これを疑似アルゴリズムで書くと以下のようになります。

matrixMultiplication()
  for i = 1 to n
    m[i][j] = 0

  for l = 2 to n
    for i = 1 to n - l + 1
      j = i + l - 1
      m[i][j] = INFTY
      for k = i to j - 1
      m[i][j] = min(m[i][j],m[i][k] + m[k + 1][j] + p[i - 1] * p[k] * p[j])

計算量の考察

上記の疑似コードでは、対象とする行列の数lを2からnまで増やしながら、その範囲を指定するiとjを決めています。
さらにその中で、iからjの範囲でkの位置を決めています。
forの3重ループになっているので、$O(n^3)$になります。

例題の解答

上記をふまえて、例題の解答をc++で書きます。

#include <stdio.h>
#include <iostream>
#include <algorithm>

using namespace std;
#define N 100

int main(void){
    int n;
    int p[N + 1],m[N + 1][N + 1] = {};
    
    cin >> n;
    for(int i = 1;i <= n;++i){
        cin >> p[i - 1] >> p[i];
    }
    for(int l = 2;l <= n;++l){
        for(int i = 1;i <= n - l + 1;++i){
            int j = i + l - 1;
            m[i][j] = __INT_MAX__;
            for(int k = i;k <= j - 1;++k){
                m[i][j] = min(m[i][j],m[i][k] + m[k + 1][j] + p[i - 1] * p[k] * p[j]);
            }
        }
    }
    cout << m[1][n] << endl;
}

参考

プログラミングコンテスト攻略のためのアルゴリズムとデータ構造

2019年10月27日日曜日

最小コストソートのアルゴリズムと例題とc++の解答

最小コストソートの考え方を理解できるようにまとめます。

最小コストソートを利用する問題

AIZU OINLINE JUDGEの最小コストソートに関する問題を例題として利用します。

問題

重さ wi(i=0,1,...,n−1) の n 個の荷物が1列に並んでいます。
これらの荷物をロボットアームを用いて並べ替えます。
1度の操作でロボットアームは荷物 i と荷物 j を持ち上げ、それらの位置を交換することができますが、
wi+wj のコストがかかります。
ロボットアームは何度でも操作することができます。

与えられた荷物の列を重さの昇順に整列するコストの総和の最小値を求めてください。

入力

1行目に整数 n が与えられます。2行目に n 個の整数 wi(i=0,1,2,...n−1) が空白区切りで与えられます。

出力

最小値を1行に出力してください。

入力例 1

5
1 5 3 4 2

出力例 1

7

最小コストソートの考え方

例として、{4,3,2,7,1,6,5}という数列の最小コストを考えます。

この数列において、各々の数値が最終的にどの位置に置かれればいいかということに絞って、 次のような入れ替えを行う。

まず先頭の4からを3番目のindexの7と交換して、

7,3,2,4,1,6,5

続いて、7を最大index6の5と交換して

5,3,2,4,1,6,7

続いて、5をindex4の1と交換して、

1,3,2,4,5,6,7

続いて、対象となる1はすでにソートの位置に置かれているので、ここで一連の流れが一旦終了することになります。

先頭から順にソートする作業を再開する、次のindex1の3に着目してindex2の2と交換する

1,2,3,4,5,6,7

2は適切な位置に置かれているので、ここで一連の流れが終了します。

つづいて、まだソート完了フラグが立っていない、index5の6に着目しますが、
6は正しい位置に置かれているので、これでソート終了ということになります。

ここで、行った処理の流れをまとめると、

①4 → 7 → 5 → 1 → 4
②3 → 2 → 3
③6 → 6

この①②③というサイクルにおいての最小コストについて考えます。

サイクルの長さが1のケース

③のようなケースでは入れ替える必要がないので、コストは0になります。

サイクルの長さが2のケース

②のようなサイクルの長さが2のケースでは、1回の交換でそれぞれの要素が最終位置に到達するので、
対象となる2つの要素の和がコストになります。
②のケースだと3 + 2 = 5がコストになる。

サイクルの長さが3以上のケース

①のケース{4 7 5 1}を考えます。
①のケースでは、以下のようにこのサイクルにおける最小コスト1を使って、それ以外の要素を移動することで、最小コストの19になります。

cost = 1 + 5 = 6
4 7 1 5
cost = 1 + 7 = 8
4 7 5 1 
cost = 1 + 4 = 5
4 1 5 7 

総コストはcost = 6 + 8 + 5 = 19となる

つまり、このケースでは、サイクルの中で最小値を使って、それ以外の要素を移動することが、
最適な方法
となる。

最小コスト移動の総コストを考える

サイクルの中の要素を$w_i$サイクル内の要素の個数をnとすると、
コストの総和 + (最小値 * (個数 - 2))と一致するので、

$\displaystyle \sum w_i + (n - 2) * min(w_i)$

となります。

最小コスト移動を行なった場合の各サイクルの総移動コスト

サイクル①のコストは0サイクル②のコストは5、
サイクル③のコストは19なので、この数列の最小総移動コストは24ということになる。

が、この考え方には反例があります。

反例を考える

数列{2,1,8,10,7,9}を例に先ほどのアルゴリズムを適応してみます。

まず、index0の先頭2とindex1の1を交換します。

1 2 8 10 7 9

交換した1は交換する必要がないので、コスト3となりここでこのサイクルは終了

続いて、8 10 7 9についてソートします。

ここではサイクルの長さが4になるので、最小コストソートを利用します。

まず、7と9を交換(コスト16)

8 10 9 7

続いて、7と10を交換(コスト17)

8 7 9 10

最後に、7と8を交換(コスト15)

7 8 9 10

総コストは、3 + (16 + 17 + 15) = 51となります。

最小コストを他のサイクルから借りてくる方法

8 → 10 → 9 → 7 → 8のサイクルについて、7と1を交換して、
8 → 10 → 9 → 1 → 8のサイクルとして、コストの総和を考えると

28 + 2 * 1に1と7の2回分の交換を加えて、このサイクルのコストは46となり、
合計しても49となるので、先ほどのアルゴリズムよりもコストが低くなります。

つまり、サイクルの外から借りてきた要素を使って、移動した方がコストが小さくなる場合があります。

最小コストをサイクルの外から借りてきた場合の計算量を考える

サイクルの外の要素を最小値をxとすると、
$2 * (min(w_i) + x$だけコストが増加するが、
$(n - 1) * (min(w_i) - x)$だけコストが減ります。
この時のコストは

$\displaystyle \sum w_i + (n - 2) * min(w_i) + 2 * (min(w_i) + x - (n - 1) * (min(w_i) - x)$ $ = \displaystyle \sum w_i + min(w_i) + (n + 1) * x$

となる。

以上を考えると、各サイクルのコストは、要素全体の最小値を借りた場合と借りない場合を計算して、
コストが小さい方を採用する必要があります。

最小コスト問題の実装のポイント

上記の仕組みを利用して、プログラムを組む上でのポイントをまとめます。

ソートされた数値がどの位置に置かれるか保存する配列を持つ

ソートされた数列の各々の数値がどの位置に入るかを記録することで、
数値を交換する動作を計算なしで行うことができます。

例えば、以下のようなコードを考えます。

#define MAX 1000
// 最大入力値
#define VMAX 10001

int n;
//数列入力用
int A[MAX];
// 入力数列の中の最小値
int s;
// 入力したAをソートした配列
int B[MAX];
// 数値をINDEXとしてソートした数列がどの位置に置かれるのかを記録する
int T[VMAX];

以下のように、ソートされた数列の各々の数値が何番目に置かれるのかを記録します。
計数ソートの考え方と同じです。

for(int i = 0;i < n;++i){
    B[i] = A[i];
}
// Bをソートする
sort(B,B + n);

for(int i = 0;i < n;++i){
    // ソートされた並び順ごとに順番を入れておく
    // ソートされた数列が1 3 5 7 9とするとT[1] = 0,T[3] = 1のように順番が入る
    T[B[i]] = i;
}

配列T[]は交換するサイクルにおいて、以下のように使えます。

for(int i = 0;i < n;++i){
    if(V[i]) continue;
    // 初回は0が対象になる
    int cursor = i;
    // 1サイクルでかかるコスト
    int S = 0;
    // 1サイクルにおける最小コスト、取りうる最大数値にしておく
    int m = VMAX;
    // サイクルの長さ
    int an = 0;
    // 1サイクル処理
    while(1){
        // ソート済みフラグ
        V[cursor] = true;
        an++;
        // 対象のコスト
        int v = A[cursor];
        // 最小コスト
        m = min(m,v);
        // 最小の値においての交換を加味しないサイクルの合計コスト
        S += v;
        // 対象が置かれるべき位置を割り出すT[v]に置かれるべき位置を入れる
        cursor = T[v];
        // サイクルが終了していたら、次に移る
        if(V[cursor]) break;
    }
    // サイクル内の最小値を利用したコスト計算(Σwi + (n - 2) * min(wi)
    int cost1 = S + (an - 2) * m;
    // 数列の最小コストを借りた場合のコスト(Σwi + min(wi) + (n + 1) * min(全数列)
    int cost2 = S + m + (an + 1) * s;
    // 最小値のコストを外から借りた場合とコストを比較する
    ans += min(cost1,cost2);
}

全体のソースコード

c++の全体のコードを載せます。

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

#define MAX 1000
// 最大入力値
#define VMAX 10001

int n;
//数列入力用
int A[MAX];
// 入力数列の中の最小値
int s;
// 入力したAをソートした配列
int B[MAX];
// 数値をINDEXとしてソートした数列がどの位置に置かれるのかを記録する
int T[VMAX];

int solve(){
    // コスト
    int ans = 0;
    // 交換終了フラグ
    bool V[MAX] = {};
    for(int i = 0;i < n;++i){
        B[i] = A[i];
    }
    // Bをソートする
    sort(B,B + n);
    
    for(int i = 0;i < n;++i){
        // ソートされた並び順ごとに順番を入れておく
        // ソートされた数列が1 3 5 7 9とするとT[1] = 0,T[3] = 1のように順番が入る
        T[B[i]] = i;
    }
    
    for(int i = 0;i < n;++i){
        if(V[i]) continue;
        // 初回は0が対象になる
        int cursor = i;
        // 1サイクルでかかるコスト
        int S = 0;
        // 1サイクルにおける最小コスト、取りうる最大数値にしておく
        int m = VMAX;
        // サイクルの長さ
        int an = 0;
        // 1サイクル処理
        while(1){
            // ソート済みフラグ
            V[cursor] = true;
            an++;
            // 対象のコスト
            int v = A[cursor];
            // 最小コスト
            m = min(m,v);
            // 最小の値においての交換を加味しないサイクルの合計コスト
            S += v;
            // 対象が置かれるべき位置を割り出すT[v]に置かれるべき位置を入れる
            cursor = T[v];
            // サイクルが終了していたら、次に移る
            if(V[cursor]) break;
        }
        // サイクル内の最小値を利用したコスト計算(Σwi + (n - 2) * min(wi)
        int cost1 = S + (an - 2) * m;
        // 数列の最小コストを借りた場合のコスト(Σwi + min(wi) + (n + 1) * min(全数列)
        int cost2 = S + m + (an + 1) * s;
        // 最小値のコストを外から借りた場合とコストを比較する
        ans += min(cost1,cost2);
    }
    return ans;
}

int main(){
    cin >> n;
    // 最小値
    s = VMAX;
    for(int i = 0;i < n;++i){
        cin >> A[i];
        s = min(s,A[i]);
    }
    int ans = solve();
    cout << ans << endl;
}

参考

プログラミングコンテスト攻略のためのアルゴリズムとデータ構造

2019年10月23日水曜日

c++で8パズルを単純な幅優先探索で解く方法

aizu online judgeから問題を引用します。

8パズルの問題

8 パズルは1つの空白を含む 3×3のマス上に 8 枚のパネルが配置され、
空白を使ってパネルを上下左右にスライドさせ、絵柄を揃えるパズルです。

この問題では、次のように空白を 0、各パネルを 1 から 8 の番号でパズルを表します。

1 3 0
4 2 5
7 8 6

1 回の操作で空白の方向に1つのパネルを移動することができ、ゴールは次のようなパネルの配置とします。

1 2 3
4 5 6
7 8 0

8パズルの初期状態が与えられるので、ゴールまでの最短手数を求めるプログラムを作成してください。

8パズルの解法

状態の遷移状態がそれほど多くないことから、単純な深さ優先探索や幅優先探索で解くことができる。
これらの解法の練習になります。

パズルの状態を管理する構造体・クラスを作る

諸々必要な変数を定義します

// 要素数
#define N 3
#define N2 N * N

// 上,右,下,左の順でチェックする
static const int dx[4] = {0,1,0,-1};
static const int dy[4] = {-1,0,1,0};
static const char dir[4] = {'u','r','d','l'};

struct Puzzle{
    // 各々のタイルの状態(spaceは0だけど、9に換算する)
    int tiles[N2];
    // 空白の位置を管理indexで0 - 8まで
    int space;
    // 解法までのルートを記録
    string path;
    
    // 同じ状態遷移か判定する
    bool operator < (const Puzzle &p) const {
        for(int i = 0;i < N2;++i){
            if(tiles[i] == p.tiles[i]) continue;
            // 最短を求めるためにタイル数値が大きいと一致していると判定する(上から1,2,3と並ぶので)
            return tiles[i] > p.tiles[i];
        }
        return false;
    }
};

ポイント

タイルの表現方法

タイルは配列で管理します。
indexを割り振ると

1[0] 2[1] 3[2]
4[3] 5[4] 6[5]
7[6] 8[7] s[8]

となります。

枝刈りの方法

mapを使用して、その状態を保存するかどうか判定します。
最短経路を求める必要があるので、保存したマップの中でindexが小さいタイルの数値が大きい場合は、
非効率とみなして、一致していると判定する
ようにします。

パズルが一致しているかチェックする関数

indexで0から保存しているので、+1してチェックしています

/**
    パズルの一致チェック
*/
bool isCorrect(Puzzle p){
    // パズルが順番通り一致しているかチェック
    for(int i = 0;i < N2;++i){
        // 1,2,3,4,5,6,7,8の順で並ぶ
        if(p.tiles[i] != (i + 1)){
          return false;
        }
    }
    return true;
}

幅優先探索を行う関数

// 幅優先探索を行う
string bfs(Puzzle p){
    // 幅優先なのでqueue
    queue Q;
    // 同じ遷移状態がないか判定するためにmapを利用
    map V;
    Puzzle u,v;
    p.path = "";
    Q.push(p);
    // 最初の状態遷移を保存
    V[p] = true;
    
    // 幅優先探索
    while(!Q.empty()){
        u = Q.front();
        Q.pop();
        
        // 正解チェック
        if(isCorrect(u)){
            return u.path;
        }
        // 空白の位置を算出する
        int sx = u.space % N;
        int sy = u.space / N;
        
        // スペースの周りの四方のタイルを動かす
        for(int r = 0;r < 4;++r){
            // 動かすタイルのindexを算出する(上,右,下,左の順)
            int tx = sx + dx[r];
            int ty = sy + dy[r];
            // 範囲外チェック
            if(tx < 0 || ty < 0 || tx >= N || ty >= N){
                continue;
            }
            v = u;
            // 入れ替え
            swap(v.tiles[u.space],v.tiles[ty * N + tx]);
            // 入れ替えた場所がspaceになる
            v.space = (ty * N) + tx;
            
            // 1.同じ状態がでなかったら
            if(!V[v]){
                V[v] = true;
                v.path += dir[r];
                Q.push(v);
            }
        }
    }
    return "unsolvable";
}

ポイント

ポイントは、1のmapを使って、最短の経路を探しているところです。
探索するか探索しないかの状態チェックを行うのに先ほどのPuzzle構造体で定義したoperator関数を使っています。

8マップのように探索量がそれほどでもないなら、これでいけるようです。

全体のソースコード

最後に全体のソースコードを載せます

#include <iostream>
#include <cmath>
#include <string>
#include <map>
#include <queue>

using namespace std;

// 要素数
#define N 3
#define N2 N * N

// 上,右,下,左の順でチェックする
static const int dx[4] = {0,1,0,-1};
static const int dy[4] = {-1,0,1,0};
static const char dir[4] = {'u','r','d','l'};

struct Puzzle{
    // 各々のタイルの状態(spaceは0だけど、9に換算する)
    int tiles[N2];
    // 空白の位置を管理indexで0 - 8まで
    int space;
    // 解法までのルートを記録
    string path;
    
    // 同じ状態遷移か判定する
    bool operator < (const Puzzle &p) const {
        for(int i = 0;i < N2;++i){
            if(tiles[i] == p.tiles[i]) continue;
            // 最短を求めるためにタイル数値が大きいと一致していると判定する(上から1,2,3と並ぶので)
            return tiles[i] > p.tiles[i];
        }
        return false;
    }
};

/**
    パズルの一致チェック
*/
bool isCorrect(Puzzle p){
    // パズルが順番通り一致しているかチェック
    for(int i = 0;i < N2;++i){
        // 1,2,3,4,5,6,7,8の順で並ぶ
        if(p.tiles[i] != (i + 1)){
          return false;
        }
    }
    return true;
}

// 幅優先探索を行う
string bfs(Puzzle p){
    // queueを使って、パズルの状態遷移を管理する
    queue Q;
    // 同じ遷移状態がないか判定するためにmapを定義する
    map V;
    Puzzle u,v;
    p.path = "";
    Q.push(p);
    // 最初の状態遷移を保存
    V[p] = true;
    
    // 幅優先探索
    while(!Q.empty()){
        u = Q.front();
        Q.pop();
        
        // 正解チェック
        if(isCorrect(u)){
            return u.path;
        }
        // 空白の位置を算出する
        int sx = u.space % N;
        int sy = u.space / N;
        
        // スペースの周りの四方のタイルを動かす
        for(int r = 0;r < 4;++r){
            // 動かすタイルのindexを算出する(上,右,下,左の順)
            int tx = sx + dx[r];
            int ty = sy + dy[r];
            // 範囲外チェック
            if(tx < 0 || ty < 0 || tx >= N || ty >= N){
                continue;
            }
            v = u;
            // 入れ替え
            swap(v.tiles[u.space],v.tiles[ty * N + tx]);
            // 入れ替えた場所がspaceになる
            v.space = (ty * N) + tx;
            
            // 同じ状態がでなかったら
            if(!V[v]){
                V[v] = true;
                v.path += dir[r];
                Q.push(v);
            }
        }
    }
    return "unsolvable";
}

int main(){
    Puzzle in;
    
    for(int i = 0;i < N2;++i){
        cin >> in.tiles[i];
        if(in.tiles[i] == 0){
            // 空白は9として管理
            in.tiles[i] = N2;
            in.space = i;
        }
    }
    
    string answer = bfs(in);
    cout << answer.size() << endl;
    
    return 0;
}

2019年9月24日火曜日

四角形内に円が内包する判定するプログラム

四角形内に縁が内包しているか判定するプログラムについてまとめたいと思います。

題材にするのは、AIZU ONINE JUDGEのCircle In Rectangleという問題です。

四角形内に円が内包するか判定する例題

長方形の中に円が含まれるかを判定するプログラムを作成してください。
次のように、長方形は左下の頂点を原点とし、右上の頂点の座標 (W,H) が与えられます。
また、円はその中心の座標 (x,y) と半径 r で与えられます。

入力

5つの整数 W、H、x、y、r が空白区切りで1行に与えられます。

出力

円が長方形の内部に含まれるなら Yes と、一部でもはみ出るならば No と1行に出力してください。

四角形内に円が内包する条件

下図より、四角形を形成する各々の直線と、円が交差(はみ出している場合)は円は内包されていません。

また、下図のように、四角形を形成する直線と円が交差してなくとも、円の中心が四角形内にない場合は、円に内包されていません。

これらの条件を組み合わせれば、四角形に円が内包されているか判定することができそうです。

円が四角形を形成する各々の線分と衝突しているか判定する方法

円の中心と直線との距離dを算出して、dと半径を比較して、半径よりもdが短かった場合に、
円と直線が衝突していると判定できます。

点と直線の距離については、こちらにまとめました。

円の中心が四角形に存在するか判定する方法

四角形が原点(0,0),width,heightのみを与えれているので、
点がwidth,heightの範囲内か判定すればokです。

点と直線の距離の公式をプログラムで表す

以下のように点と直線の距離を返す関数を定義しました。

struct Vec2{
    int x;
    int y;
    
    Vec2(int x1,int y1){
        x = x1;
        y = y1;
    }
    Vec2(){
        
    }
};

float getlengthLineBetweenPoint(float a,float b,float c,Vec2 point){
    // 点と直線の距離の公式に当てはめる
    return abs((a * point.x + b * point.y + c) / sqrt((a * a) + (b * b)));
}

全体のコード

問題の答えとなるコードは以下のように書きました。

#include <iostream>
#include <stdlib.h>
#include <math.h>

using namespace std;

struct Vec2{
    int x;
    int y;
    
    Vec2(int x1,int y1){
        x = x1;
        y = y1;
    }
    Vec2(){
        
    }
};

struct Line{
    // vector2D
    Vec2 base;
    // vector2D
    Vec2 direction;
};

float getlengthLineBetweenPoint(float a,float b,float c,Vec2 point){
    // 点と直線の距離の公式に当てはめる
    return abs((a * point.x + b * point.y + c) / sqrt((a * a) + (b * b)));
}

int main(void){
    int width,height,centerX,centerY,radius;
    cin >> width >> height >> centerX >> centerY >> radius;
    
    // 原点を通る縦直線(x = 0)
    float distance1 = getlengthLineBetweenPoint(1,0,0,Vec2(centerX,centerY));
    // 原点を通る横直線(y = 0)
    float distance2 = getlengthLineBetweenPoint(0,1,0,Vec2(centerX,centerY));
    // 縦直線
    float distance3 = getlengthLineBetweenPoint(1,0,-width,Vec2(centerX,centerY));
    // 横直線
    float distance4 = getlengthLineBetweenPoint(0,1,-height,Vec2(centerX,centerY));
    
    string str = "No";
    // 円の半径の距離が全ての直線との距離未満ならばぶつかっていないことになる
    // 円と直線がぶつかっていない
    if(distance1 >= radius && distance2 >= radius && distance3 >= radius && distance4 >= radius){
        // かつ円の中心が四角形の内部になる
        // x座標
        if(centerX > 0 && centerX < width){
            // y座標
            if(centerY > 0 && centerY < height){
                str = "Yes";
            }
        }
    }
    cout << str << endl;
    
    return 0;
}

まとめ

四角形に円が内包するかの条件は、
1.円と各々の直線が交差しているか
2.四角形内に円の中心が内包しているか

で判定することができる。

効率の良いやり方は、いくらでもあると思うけど。
単純でできそうなやり方をいくつか考えたけれど、一般的な方法で落ち着いてしまった。。。

2019年8月27日火曜日

隣接リストを使用した幅優先探索でALDS1_11_Cを解く

以前は隣接行列と幅優先探索を使って、ALDS_1_11_Cの問題を解きましたが、 今回は、隣接リストと幅優先探索を使って、この問題を解いてみたいと思います。

隣接リストを使用した幅優先探索・深さ優先探索の計算量は、こちらに書いたとおり、
O(|V| + |E|)の計算量になるので、隣接行列を使用した時よりも計算量が少なくなります。
(V = 頂点数,E = Vに対して隣接する頂点)

隣接リストを表現する方法

隣接リストを表現するために、c++のvectorを使います。

隣接リストの入力は以下のような形式で渡されるので、

1 2 2 4

以下のように隣接リストを定義して、挿入すればよいです。

// 隣接リストの定義
vectorlist[N];

    // nは頂点数
    for(int i = 0;i < n;++i){
        // uは頂点番号
        // kはuに隣接する頂点の数
        cin >> u >> k;
        // 頂点を移動する(便宜上0から始まるようにずらす)
        u--;
        // 頂点の数だけ回す
        for(int j = 0;j < k;++j){
            cin >> v;
            // 隣接行列に挿入する
            list[u].push_back(v - 1);
        }
    }

隣接リストによる幅優先探索のやり方

では、隣接リストにより幅優先探索を行う箇所を見ます。

隣接行列では、行列の数列分探索していましたが、隣接リストでは、頂点が隣接している分だけ回すだけでokです。

    while(!q.empty()){
        u = q.front();
        q.pop();
        // 頂点分探索する
        for(int i = 0;i < list[u].size();++i){
            // 隣接する頂点
            int v = list[u][i];
            // 頂点が未訪問
            if(d[v] == INFTY){
                // 対象の頂点の距離は前の頂点 + 1
                d[v] = d[u] + 1;
                q.push(v);
            }
        }
    }

全体のソースコード

最後に全体のソースコードを載せます。

#include<iostream>
#include<queue>
#include<vector>
#include<stack>

using namespace std;

#define N 100
#define INFTY 99999999
vectorlist[N];

int n;

// 距離 colorも同時に管理する
int d[N];

// 幅優先探索
void bfs(int s){
    queue q;
    q.push(s);

    // 距離を初期化
    for(int i = 0;i < n;++i){
        d[i] = INFTY;
    }
    d[s] = 0;
    // 対象の頂点番号
    int u;

    while(!q.empty()){
        u = q.front();
        q.pop();
        // 頂点分探索する
        for(int i = 0;i < list[u].size();++i){
            // 隣接する頂点
            int v = list[u][i];
            // 頂点が未訪問
            if(d[v] == INFTY){
                // 対象の頂点の距離は前の頂点 + 1
                d[v] = d[u] + 1;
                q.push(v);
            }
        }
    }
}

int main(){
    int u,k,v;
    
    cin >> n;
    
    for(int i = 0;i < n;++i){
        cin >> u >> k;
        // 頂点を移動する
        u--;
        // 頂点の数だけ回す
        for(int j = 0;j < k;++j){
            cin >> v;
            // 隣接行列に挿入する
            list[u].push_back(v - 1);
        }
    }
    
    bfs(0);
    
    for(int i = 0;i < n;++i){
        cout << i + 1 << " " << ((d[i] == INFTY) ? (-1) : d[i]) << endl;
    }
    
    return 0;
}

2019年8月26日月曜日

グラフにおける隣接リストの表現方法

グラフにおける隣接リストの表現方法についてまとめたいと思います。
まずは、隣接リストの定義についてまとめます。(AIZU ONLINE JUDGEから引用)

隣接リストとは

隣接リスト表現では、V の各頂点に対して1個、合計 |V| 個のリスト Adj[|V|] でグラフを表します。
頂点 u に対して、隣接リスト Adj[u] は E に属する辺 (u,vi) におけるすべての頂点 vi を含んでいます。
つまり、Adj[u] は G において u と隣接するすべての頂点からなります。

この表現方法だと一目で理解するのは、難しいですが、
隣接リストの入力方法を見れば、理解しやすいです。
隣接リストの入力方法では、以下のような形で与えられます。

隣接リストの入力方法

最初の行に G の頂点数 n が与えられます。
続く n 行で各頂点 u の隣接リスト Adj[u] が以下の形式で与えられます:

u k v1 v2 ... vk

u は頂点の番号、k は u の出次数、v1v2...vk は u に隣接する頂点の番号を示します。

つまり、隣接リストとして以下のような形で入力されます。

1 2 2 4

これは頂点1に2個隣接している頂点があり、それが2と4であることを示しています。

この形式を表現するには、vectorが適しています。

隣接リストをvectorで表現する

以下のようにvectorの領域を頂点分用意することによって、
隣接行列を簡単に表現することができます。

// 頂点分vectorの領域を用意する
vectorgraph[n];
// 頂点uに隣接するvを追加する
graph[u].push_back(v);

// 頂点uに隣接する頂点を探索する場合
for(int i = 0;i < graph[u].size;++i){
 // 隣接する頂点を取得する
 int v = graph[u][i];
} 

隣接リストを用いた計算量

隣接リストを用いて、幅優先探索・深さ優先探索は、
各頂点を一度訪問して、隣接リストの頂点を調べるので、
O(|V| + |E|)の計算量になります。

2019年8月23日金曜日

幅優先探索を使用してグラフの頂点から頂点への最短距離を求める方法をALDS_1_11_Cを使って学ぶ

キューを使った幅優先探索をAIZU ONLINE JUDGEの例題を利用して、学びたいと思います。
幅優先探索の概念についてはこちらにまとめています。

では、例題を見ていきます。

例題

与えられた有向グラフ G=(V,E) について、
頂点 1 から各頂点への最短距離 d(パス上の辺の数の最小値)を求めるプログラムを作成してください。
各頂点には 1 から n までの番号がふられているものとします。
頂点 1 からたどり着けない頂点については、距離として-1 を出力してください。

入力

最初の行に G の頂点数 n が与えられます。続く n 行で各頂点 u
の隣接リストが以下の形式で与えられます:

u k v1 v2 ... vk

uは頂点の番号、k は u の出次数、
v1v2...vk は u に隣接する頂点の番号を示します。

制約

1≤n≤100

出力

各頂点について id、d を1行に出力してください。
id は頂点の番号、d は頂点 1 からその頂点までの距離を示します。頂点番号順に出力してください。

解法

幅優先探索は、始点sから距離k + 1にある頂点を探索する前に、距離kの頂点をすべて発見しているので、
始点から各頂点までの最短距離を順番に求めることができます。

幅優先探索のアルゴリズム

1.始点sをキューに入れる

キューが空でない限り以下の操作を繰り返す

  • キューから頂点uを取り出して訪問する(訪問完了)
  • 頂点uに隣接し、未訪問の頂点vについて、始点から頂点までの距離を表す変数d[v]をd[u] + 1と更新し、
    頂点vをキューに入れる

グラフを例に幅優先探索が行われる流れを見る

下図のグラフを例として、頂点0から幅優先探索が行われる流れを追います。

1.始点である頂点0をキューに挿入します。
始点0からの距離を0とします。

0

2.キューから0を取り出して、訪問します。
頂点0に隣接して未訪問の1,3,5の頂点をキューに挿入します。 頂点1,3,5までの距離は頂点0までの距離+1となります。

1,3,5

3.キューの先頭から1を取り出して、訪問します。
頂点1に隣接する未訪問の頂点はないので、次に進みます。

3,5

4.キューの先頭から3を取り出して、訪問します。
頂点3に隣接する未訪問の頂点は、2のみなので、2をキューに挿入します。 頂点2までの距離は、頂点3までの距離+1になります。

5,2

5.キューの先頭から5を取り出して、訪問します。
頂点5に隣接する未訪問の頂点4,6をキューに入れます。
頂点4,6までの距離は頂点5までの距離+1になります。

2,4,6

6.キューの先頭から2を取り出して、訪問します。
頂点2に隣接する未訪問の頂点はないので、次に進みます。

4,6

7.キューの先頭から4を取り出して、訪問します。
頂点4に隣接する未訪問の頂点はないので、次に進みます。

6

8.キューの先頭から6を取り出して、訪問します。
頂点6に隣接する未訪問の頂点はないので、次に進みます。

キューが空になったので、幅優先探索終了です。
探索終了時に頂点0から各々の頂点への最短距離が求められます。

幅優先探索に必要な変数

キューを利用した深さ優先探索に必要な変数は以下のようになります。

変数名例 役割
color[N] 未訪問・訪問中・訪問済みの状態を表すための変数
M[n][n] グラフの隣接行列を表すための変数、頂点が隣接している場合は1(true)など
Queue q 使用するキュー
d[n] 始点からインデックスを頂点番号とする距離を保存しておく

隣接行列を利用した幅優先探索の疑似アルゴリズム

bfs()

 各頂点についてcolorを未訪問に設定する
 各頂点において、d[u]を到達できないことを示す値を入れる
 // 始点をsとします
 color[s] = 訪問中
 d[s] = 0
 q.enque(s)
 
 while Qが空でない
  u = Q.dequeue()
  for vが0から|V| - 1まで
   if M[u][v] && color[v] == 未訪問
    color[v] = 訪問中
    d[v] = d[u] + 1
    Q.enqueue(v)
  color[u] = 訪問済み
  
 

計算量考察

隣接行列を用いた幅優先探索は、各頂点について、すべての頂点に隣接しているかどうかを調べるので、
$O(|V|^2)$のアルゴリズムになります。

c++を用いた解答例

では最後に隣接行列で幅優先探索を行う解答例を載せます


#include<iostream>
#include<queue>

using namespace std;

#define N 100
#define INFTY 99999999

int n,
// 隣接行列
M[N][N];
// 距離 colorも同時に管理する
int d[N];

void bfs(int s){
    queue q;
    q.push(s);
    
    // 距離を初期化
    for(int i = 0;i < n;++i){
        d[i] = INFTY;
    }
    d[s] = 0;
    // 対象の頂点番号
    int u;
    
    while(!q.empty()){
        u = q.front();
        q.pop();
        
        for(int v = 0;v < n;++v){
            // 頂点が隣接しているかつ未訪問
            if(M[u][v] == 1 && d[v] == INFTY){
                // 対象の頂点の距離は前の頂点 + 1
                d[v] = d[u] + 1;
                q.push(v);
            }
        }
    }
    
}

int main(){
    int u,k,v;
    
    cin >> n;
    
    // 隣接行列を初期化
    for(int i = 0;i < n;++i){
        for(int j = 0;j < n;++j){
            M[i][j] = 0;
        }
    }
    
    for(int i = 0;i < n;++i){
        cin >> u >> k;
        // 頂点を移動する
        u--;
        // 頂点の数だけ回す
        for(int j = 0;j < k;++j){
            cin >> v;
            --v;    // ずらす
            M[u][v] = 1;
        }
    }
    
    bfs(0);

    for(int i = 0;i < n;++i){
        cout << i + 1 << " " << ((d[i] == INFTY) ? (-1) : d[i]) << endl;
    }

    return 0;
}

2019年8月22日木曜日

スタックによる深さ優先探索をALDS1_11_Bを使って学ぶ

スタックを使用して、深さ優先探索を行う方法をAIDU_ONLINE_JUDGEのALDS1_11_Bを利用して学びたいと思います。

例題

深さ優先探索(Depth First Search: DFS)は、可能な限り隣接する頂点を訪問する、という戦略に基づくグラフの探索アルゴリズムです。
未探索の接続辺が残されている頂点の中で最後に発見した頂点 v の接続辺を再帰的に探索します。

v の辺をすべて探索し終えると、vを発見したときに通ってきた辺を後戻りして探索を続行します。
探索は元の始点から到達可能なすべての頂点を発見するまで続き、 未発見の頂点が残っていれば、その中の番号が一番小さい1つを新たな始点として探索を続けます。

深さ優先探索では、各頂点に以下のタイムスタンプを押します:

  • タイムスタンプ d[v]: vを最初に発見した発見時刻を記録します。
  • タイムスタンプ f[v]: vの隣接リストを調べ終えた完了時刻を記録します。

以下の仕様に基づき、与えられた有向グラフ G=(V,E)
に対する深さ優先探索の動作を示すプログラムを作成してください:

Gは隣接リスト表現の形式で与えられます。
各頂点には 1 から nまでの番号がふられています。 各隣接リストの頂点は番号が小さい順に並べられています。 プログラムは各頂点の発見時刻と完了時刻を報告します。 深さ優先探索の過程において、訪問する頂点の候補が複数ある場合は頂点番号が小さいものから選択します。 最初に訪問する頂点の開始時刻を 1 とします。

入力

最初の行に G の頂点数 n が与えられます。続く n 行で各頂点 u の隣接リストが以下の形式で与えられます:

u k v1 v2 ... vk

u は頂点の番号、k は u の出次数、v1v2...vk  は u に隣接する頂点の番号を示します。

出力

各頂点について id、 d、 fを空白区切りで1行に出力してください。
id は頂点の番号、d はその頂点の発見時刻、f はその頂点の完了時刻です。頂点の番号順で出力してください。

制限

1≤n≤100

入力例

4
1 1 2
2 1 4
3 0
4 1 3

出力例

1 1 8
2 2 7
3 4 5
4 3 6

スタックを利用した深さ優先探索

スタックにまだ探索中の頂点を一時的に保存しておくことにより、
深さ優先探索を行います。

手順は以下のようになります。

1.一番最初に訪問する頂点をスタックに入れる

2.スタックに頂点が積まれている限り、以下の処理を繰り返します。

  • スタックの先頭に積まれている頂点uを訪問する
  • 現在訪問中の頂点uから次の頂点vへ探索を移行するときに、
    頂点vをスタックに積みます。
    現在訪問中の頂点uに未訪問の隣接する頂点なければuをスタックから解放します。

スタックの処理例

以下のグラフを例にスタックによる、深さ優先探索を行って見たいと思います。

1.最初に訪問する頂点を0としてスタックに積みます

2.スタックの先頭の0を訪問します。
0に隣接し、かつ未訪問で、数の小さい方の頂点1をスタックに積みます。

3.スタックの先頭1を訪問します。
1に隣接しており、未訪問の3をスタックに積みます。

4.スタックの先頭3を訪問します。
3に隣接し、未訪問の2をスタックに積みます

5.スタックの先頭2を訪問します。
頂点2に隣接して未訪問な頂点は存在しないので、2をスタックから削除します。

6.スタック先頭の3を再び訪問します。 頂点3に隣接して、未訪問の頂点6をスタックに積みます。

7.スタックの先頭6を訪問します。
頂点6に隣接して、未訪問の頂点5をスタックに積みます。

8.スタックの先頭5を訪問します。
頂点5に隣接して、未訪問の頂点4をスタックに積みます。

8.スタックの先頭4を訪問します。
頂点4に隣接して、未訪問の頂点の頂点はないので、4をスタックから削除します。

9.スタックの先頭5訪問します。
頂点5に隣接して、未訪問の頂点の頂点はないので、5をスタックから削除します。

10.スタックの先頭6訪問します。
頂点6に隣接して、未訪問の頂点の頂点はないので、6をスタックから削除します。

11.スタックの先頭3訪問します。
頂点3に隣接して、未訪問の頂点の頂点はないので、3をスタックから削除します。

12.スタックの先頭1訪問します。
頂点1に隣接して、未訪問の頂点の頂点はないので、1をスタックから削除します。

13.スタックの先頭0訪問します。
頂点0に隣接して、未訪問の頂点の頂点はないので、0をスタックから削除します。

14.スタックが空となり、すべての頂点の訪問が完了しました。

スタックによる深さ優先探索に必要な変数一覧

スタックによる深さ優先探索に必要な変数は以下のようになります。

変数名例 役割
color[N] 未訪問・訪問中・訪問済みの状態を表すための変数
M[n][n] グラフの隣接行列を表すための変数、頂点が隣接している場合は1(true)など
Stack s 訪問中の頂点を退避しておくためのスタック

スタックによる深さ優先探索の疑似アルゴリズム

以上の処理を疑似アルゴリズムにまとめてみます。

dfs(u)
 S.push(u) // 始点をスタックに追加
 color[u] = 訪問中
 d[u] = ++time; // 訪問時間を代入する
 
 while(sが空でない)
  u = S.top()
  v = next(u) // uに隣接している頂点を昇順になるように取得する
  if v != NIL
   if color[v] == 未訪問
    color[v] = 訪問中
    d[v] = ++time; // 訪問時間を記録
    S.push(v)
   // 隣接する頂点がない場合   
  else
   S.pop()
   color[v] = 訪問済み
   f[u] = ++time; // 探索が終わった時間を記録

stackを使い隣接行列を用いた計算量

各々の頂点について、隣接しているすべての頂点を調べるので、
計算量は$O[n^2]$の計算量になります。

c++による、stackを利用した例題の解法

最後に隣接行列とstackを利用した問題のc++での解法を載せます。


#include<iostream>
#include<stack>
using namespace std;

#define N 100
// 未訪問
#define WHITE 0
// 訪問中
#define GRAY 1
// 訪問済み
#define BLACK 2

int n,M[N][N];
int color[N];
// 訪問時間
int d[N];
// 探索終了時間
int f[N];
// 時間
int t = 0;
// 探索済みの頂点を一から探索しないように、nt[u]にすでに探索済みの番号を入れ、
// 次回は探索した移行の頂点から探索されるようにする
int nt[N];

/**
    uに隣接するvを番号順に取得する
*/
int next(int u){
    // 探索済み移行の頂点から探索を始める
    for(int v = nt[u];v < n;++v){
        nt[u] = v + 1;
        if(M[u][v]){
            return v;
        }
    }
    return -1;
}

/**
    スタックを用いた深さ優先探索
*/
void dfs(int u){
    for(int i = 0;i < n;++i){
        nt[i] = 0;
    }
    
    stack S;
    S.push(u);
    color[u] = GRAY;    // 訪問中
    d[u] = ++t;      // 訪問時間
    
    while(!S.empty()){
        int topU = S.top();
        // 隣接する頂点
        int v = next(topU);
        // 隣接する頂点が見つかった場合
        if(v != -1){
            // 未探索のケース
            if(color[v] == WHITE){
                color[v] = GRAY;
                // 訪問時間
                d[v] = ++t;
                S.push(v);
            }
        }else{
            // 未訪問の頂点が見つからなかった
            S.pop();
            // 探索終了
            color[topU] = BLACK;
            // 探索終了時間
            f[topU] = ++t;
        }
    }
}

int main(){
    int u,k,v;
    cin >> n;
    
    // 初期化
    for(int i = 0;i < n;++i){
        color[i] = WHITE;
        for(int j = 0;j < n;++j){
            M[i][j] = 0;
        }
    }
    
    for(int i = 0;i < n;++i){
        cin >> u >> k;
        // 頂点を移動する
        u--;
        // 頂点の数だけ回す
        for(int j = 0;j < k;++j){
            cin >> v;
            --v;    // ずらす
            M[u][v] = 1;
        }
    }
    for(int i = 0;i < n;++i){
        if(color[i] == WHITE){
            dfs(i);
        }
    }
    
    for(int i = 0;i < n;++i){
        cout << i + 1 << " " << d[i] << " " << f[i] << endl;
    }
    
    return 0;
}

2019年8月21日水曜日

再帰関数による深さ優先探索をALDS1_11_Bを例題にして学ぶ

再帰関数による深さ優先探索をALDS1_11_Bを例題にして、学びます。
深さ優先探索の考え方については、こちらにまとめています。

例題

深さ優先探索(Depth First Search: DFS)は、可能な限り隣接する頂点を訪問する、という戦略に基づくグラフの探索アルゴリズムです。
未探索の接続辺が残されている頂点の中で最後に発見した頂点 v の接続辺を再帰的に探索します。

v の辺をすべて探索し終えると、vを発見したときに通ってきた辺を後戻りして探索を続行します。
探索は元の始点から到達可能なすべての頂点を発見するまで続き、 未発見の頂点が残っていれば、その中の番号が一番小さい1つを新たな始点として探索を続けます。

深さ優先探索では、各頂点に以下のタイムスタンプを押します:

  • タイムスタンプ d[v]: vを最初に発見した発見時刻を記録します。
  • タイムスタンプ f[v]: vの隣接リストを調べ終えた完了時刻を記録します。

以下の仕様に基づき、与えられた有向グラフ G=(V,E)
に対する深さ優先探索の動作を示すプログラムを作成してください:

Gは隣接リスト表現の形式で与えられます。
各頂点には 1 から nまでの番号がふられています。 各隣接リストの頂点は番号が小さい順に並べられています。 プログラムは各頂点の発見時刻と完了時刻を報告します。 深さ優先探索の過程において、訪問する頂点の候補が複数ある場合は頂点番号が小さいものから選択します。 最初に訪問する頂点の開始時刻を 1 とします。

入力

最初の行に G の頂点数 n が与えられます。続く n 行で各頂点 u の隣接リストが以下の形式で与えられます:

u k v1 v2 ... vk

u は頂点の番号、k は u の出次数、v1v2...vk  は u に隣接する頂点の番号を示します。

出力

各頂点について id、 d、 fを空白区切りで1行に出力してください。
id は頂点の番号、d はその頂点の発見時刻、f はその頂点の完了時刻です。頂点の番号順で出力してください。

制限

1≤n≤100

入力例

4
1 1 2
2 1 4
3 0
4 1 3

出力例

1 1 8
2 2 7
3 4 5
4 3 6

解法

隣接リストで与えられる入力に対して、隣接行列を作り、それをGraphとすると Graphに対し、深さ優先探索を再帰関数で行うプログラムは以下のようになります。

void visit(int i){
    // 1
    visited[i] = true;
    // 2
    for(int j = 0;j < max;++j){
        if(Graph[i][j] == 1 && !visited[j]){
            visit(j);
        }
    }
}

1.頂点を探索したことを示すフラグvisitedを用意し、
探索した段階で、フラグを立てます。

2.与えられた頂点から隣接行列を調べていき、隣接する頂点を調べます。
すでに探索している頂点の場合処理を飛ばすことで、ループになることを防いでいます。

タイムスタンプへの対応方法

問題では、頂点の探索開始時刻と探索終了時刻を記録する必要があります。

探索開始時点は、探索する前である、関数visitの2の処理の前におけばよく、
探索終了時点は、2の派生する点の探索が終わった時点で換算することになります。

c++での解法例

それでは、以上を踏まえたc++での解法例を書きます。

#include <iostream>

using namespace std;

#define N 100

int Graph[N][N];
// 頂点を発見した時刻
int d[N];
// 頂点の探索を終えた時刻
int f[N];
// かかった時間
int dTime = 0;
// 頂点を発見したかのフラグ
bool visited[N];

void visit(int i,int max){
    
    visited[i] = true;
    d[i] = ++dTime;
    for(int j = 0;j < max;++j){
        if(Graph[i][j] == 1 && !visited[j]){
            visit(j,max);
        }
    }
    // すべての頂点の探索が終了した
    f[i] = ++dTime;
}

int main(){
    int n,u,k,v;

    cin >> n;
    for(int i = 0;i < n;++i){
        cin >> u >> k;
        // 頂点を移動する
        u--;
        // 頂点の数だけ回す
        for(int j = 0;j < k;++j){
            cin >> v;
            --v;    // ずらす
            Graph[u][v] = 1;
        }
    }
    for(int i = 0;i < n;++i){
        if(!visited[i]){
            visit(i, n);
        }
    }
    
    for(int i = 0;i < n;++i){
        cout << i + 1 << " " << d[i] << " " << f[i] << endl;
    }
    
    return 0;
}

計算量を考える

隣接行列を使い2重ループですべてのパターンを調べているので、
$O(N^2)$の計算量がかかることになります。

2019年8月20日火曜日

ALDS1_11_A問題 隣接リスト・隣接行列についての解法

ALDS1_11_A問題(難易度1)の解法についてまとめます。

問題

グラフ G=(V,E)
の表現方法には隣接リスト(adjacency list) による表現と隣接行列(adjacency matrices)による表現があります。

隣接リスト表現では、Vの各頂点に対して1個、
合計 |V| 個のリスト Adj[|V|] でグラフを表します。
頂点 u に対して、隣接リスト Adj[u] は E に属する辺 (u,vi) におけるすべての頂点 vi を含んでいます。
つまり、Adj[u] は G において uと隣接するすべての頂点からなります。

一方、隣接行列表現では、頂点 i
から頂点 j へ辺がある場合 aij が 1、ない場合 0 であるような |V|×|V| の行列 A
でグラフを表します。

隣接リスト表現の形式で与えられた有向グラフ G の隣接行列を出力するプログラムを作成してください。G は n(=|V|) 個の頂点を含み、それぞれ 1 から n までの番号がふられているものとします。

入力

最初の行に G の頂点数 n が与えられます。
続く n 行で各頂点 u の隣接リスト Adj[u] が以下の形式で与えられます:

u k v1 v2 ... vk

uは頂点の番号、k は u の出次数、v1v2...vk は u に隣接する頂点の番号を示します。

出力

出力例に従い、G の隣接行列を出力してください。aij の間に1つの空白を入れてください。

制限

1≤n≤100

入力例

4
1 2 2 4
2 1 4
3 0
4 1 3

出力例

0 1 0 1
0 0 0 1
0 0 0 0
0 0 1 0

問題解説

隣接リストとして与えられた、数値を隣接行列に直して出力する問題です。
隣接リストや隣接行列を知らなくとも、解説にどういうものかが書いてあるので、
知らなくとも解くことができます。

制約にNについて100以下であると定義されているので、動的に領域を確保しなくとも、
N分の配列を用意し、隣接行列として利用し、そのまま出力に利用することができます。

隣接行列についてはこちらにまとめました。

解法のコツ

動的に要素を受け取る方法

入力を受け取る際に、

k v1 v2 ... vk

とkの数に従い動的に入力を受け取る箇所があります。
これをc++で書くと以下のようになります。

    int n,u,k,v;
  
    cin >> n;
    // nの個数分回す
    for(int i = 0;i < n;++i){
        // 一旦Kを受け取る
        cin >> u >> k;
        u--;
        // 頂点の数だけ回す
        for(int j = 0;j < k;++j){
            cin >> v;
        }
    }

あとは、問題に書いてあることをそのまま出力すればokです。

c++での解法

#include

using namespace std;

#define N 100

int M[N][N];

int main(){
    int n,u,k,v;
  
    cin >> n;
    for(int i = 0;i < n;++i){
        cin >> u >> k;
        // 頂点を移動する
        u--;
        // 頂点の数だけ回す
        for(int j = 0;j < k;++j){
            cin >> v;
            --v;    // ずらす
            M[u][v] = 1;
        }
    }
    
    for(int i = 0;i < n;++i){
        for(int j = 0;j < n;++j){
            // 2個目からは空白が必要
            if(j){
                cout << " ";
            }
            cout << M[i][j];
        }
        cout << endl;
    }
    
    return 0;
}

2019年8月18日日曜日

グラフおける隣接行列について

グラフにおける隣接行列について、概念をまとめたいと思います。

隣接行列とは

隣接行列は、2次元の配列でグラフを表現する方法です。
配列のインデックスが各々の頂点番号を表します。

隣接行列の2次元配列をMとするとM[i][j]が、頂点iと頂点jの関係を示します。
頂点と頂点が辺でつながれている場合は1やtrueで表し、
頂点と頂点間で繋がりがない場合は、0やfalseで表現します。

無効グラフにおける隣接行列

隣接行列の2次元配列をMとします。
無効グラフの隣接行列では、頂点i,jの間に辺があるとすると、
、M[i][j]、M[j][i]の両方の値を1にし、辺がない場合は0にします。 隣接行列は、右上と左下が対象になります。

下図のような無効グラフがあるとします。

これを隣接行列に直すと以下のようになります。
頂点と頂点が辺でつながれている場合は配列の要素が1になります。

有効グラフにおける隣接行列

有効グラフにおける隣接行列では、頂点iから頂点jに向かって辺がある場合、
M[i][j]の値を1とし、辺がない場合0とします。
無効グラフと違い頂点同士必ず相互に1とはなりません。

下図のような有効グラフがあるとします。

この有効グラフの隣接行列は以下のようになります。

隣接行列の長所

配列M[i][j]により、辺(u,v)を参照できるので、頂点uと頂点vの関係をO(1)で参照することができます。

隣接行列の短所

  • 一つの隣接行列では頂点uから頂点vへの関係を一つしか記録できない点
  • 頂点の数の2乗のメモリを使用する点

参考
プログラミングコンテスト攻略のためのアルゴリズムとデータ構造

2019年8月15日木曜日

AtCorderのABC087B - Coins問題を全探索で解く

AtCorderの問題ABC087B-Coinsを全探索で解く解法を書きます。
問題は以下になります。

問題

あなたは、500 円玉を A 枚、100 円玉を B 枚、50 円玉を C 枚持っています。
これらの硬貨の中から何枚かを選び、合計金額をちょうど X円にする方法は何通りありますか。

同じ種類の硬貨どうしは区別できません。2 通りの硬貨の選び方は、ある種類の硬貨についてその硬貨を選ぶ枚数が異なるとき区別されます。

制約

0≤A,B,C≤50
A+B+C≥1
50≤X≤20,000
A,B,Cは整数である
Xは 50 の倍数である

入力

入力は以下の形式で標準入力から与えられる。

A
B
C
X

出力

硬貨を選ぶ方法の個数を出力せよ。

問題解説

同じ硬貨の枚数が重複しない、組み合わせの総数を求めています。
制限を見る限り一つの硬貨につき50枚までとリミットが少ないので、
全探索で解いてみることにしました。

解法のコツ

数値の一番大きな硬貨500円を選択するかしないかの組み合わせから試すことで、
合計金額がマイナスになった時に、後半の組み合わせを無視することができるので、
大 ⇨ 小という風にコインの組み合わせをループを組んでいきます。

500円を選択するかしないかの組み合わせをチェックする場合は、
forで500円の個数枚分ループを回すとして、
最初の0は500円を未選択の状態とし、1の場合を1枚選択した状態、
2の場合を2枚選択した状態とし、目標の合計金額から枚数分減らしていきます。

合計金額から500円を引いた時に、合計金額がマイナスとなった時は、それ以降計算しても合計金額が0になることはないので、
組み合わせのループを抜けて、組み合わせの数を表示させるようにしました。

c++での解法例

#include <iostream>

using namespace std;

int main(){
    int a,b,c,sum;
    int success = 0;
    cin >> a;
    cin >> b;
    cin >> c;
    cin >> sum;
    
    // 500
    for(int i = 0;i <= a;++i){
        // 初回は選ばない
        if(i != 0){
            // 2回目以降は選択するつまり、合計金額を減らす
            sum -= 500;
            if(sum == 0){
                success++;
                i = a;
            }
            // -になった場合これ以降の組み合わせは無駄になる
            if(sum < 0){
                break;
            }
        }
        // 100
        for(int j = 0;j <= b;++j){
            // 初回は選ばない
            if(j != 0){
                sum -= 100;
                if(sum == 0){
                    success++;
                    sum += j * 100;
                    j = b;
                    break;
                }
                // 金額がマイナスになったら、合計金元に戻して500円の組み合わせから
                else if(sum < 0){
                    sum += j * 100;
                    break;
                }
            }
            // 50
            for(int k = 0;k < c;++k){
                sum -= 50;
                if(sum == 0){
                    success++;
                    sum += (k + 1) * 50;
                    k = c;
                    break;
                }
                else if(sum < 0){
                    break;
                }
                // ループが終了したら金額を戻す
                if(k + 1 >= c){
                    sum += (k + 1) * 50;
                }
            }
            // ループが終了したら金額を戻す
            if(j + 1 > b){
                sum += j * 100;
            }
        }
    }
    
    cout << success << endl;
    return 0;
}

2019年8月7日水曜日

8クイーン問題をバックトラック方法を用いて解く方法

AIDU ONLINE JUDGEの8クイーン問題を解いてみたいと思います。

8クイーンの例題

8 クイーン問題とは、8×8 のマスから成るチェス盤に、どのクイーンも他のクイーンを襲撃できないように、8 つのクイーンを配置する問題です。
チェスでは、クイーンは次のように8方向のマスにいるコマを襲撃することができます。
すでにクイーンが配置されている k 個のマスが指定されるので、
それらを合わせて8つのクイーンを配置したチェス盤を出力するプログラムを作成してください。

入力

1行目に整数 k が与えられます。続く k 行にクイーンが配置されているマスが2つの整数 r, c で与えられます。
r, c はそれぞれ0から始まるチェス盤の行と列の番号を表します。

出力

出力は 8×8 のチェス盤を表す文字列で、クイーンが配置されたマスを 'Q'、配置されていないマスを '.' で表します。

単純な計算量を考える

配置する8つのクイーンのすべての組み合わせを考えると、
8×8のマス目から8個を選ぶ組み合わせなので、
${}_64 \mathrm{ C }_8$ = 4426165368通りになります。

ここから、1つの行に1つしかクイーンを配置することができないので、$8^8 = 16777216$通りに、
さらに、2つ以上のクイーンが同じ列に配置できないことを考えると、8! = 40320通りになります。

バックトラックを用いた解法

バックトラックという方法を用いて、問題を解いてみます。
以下のように処理をします。

1.1行目の任意のマスにクイーンを配置する

2.1行目に配置したクイーンに襲撃されない位置に2行目のマスにクイーンを配置する

処理をまとめますと、

各クイーンが他のどのクイーンも襲撃しないように、i個のクイーンを最初のi行に配置した状態で、
これらのどのクイーンにも襲撃されない(i + 1)行目のマスに、クイーンを配置する

もし、そのようなマスが(i + 1)行目にないならば、i行目に戻り、i行目で行なっていた
他のクイーンに襲撃されないマス目探しを続けます。
そのようなマス目がなければ、さらに(i - 1)行に戻り、処理を続けていきます。

このように、可能性のある状態を調べていき、現在の状態からは得られないとわかった時点で
探索を打ち切り、一つ前の状態に戻って、途中から探索を再開する手法をバックトラックというそうです。

8クイーン問題のバックトラック手法に用いる変数と役割

クイーンの配置状態を記録するために、以下の変数を定義します。

変数名例 役割
row[N] rox[i]に行が襲撃されているかどうかを記録する(rowは行の意味)
col[N] col[j]に列が襲撃されているかどうかを記録する(colは列の意味)
dpos[2N - 1] dpos[i + j]に左下方向の列が襲撃されているかどうかを記録する
dneg[2N - 1] dneg[i - j + (N - 1)]にに右下方向の列が襲撃されているかどうかを記録する

下図が変数dposがどういう状態を表すかを示した図になります。
斜めに領域を取るので、dposの容量は2N - 1になります。
この図を左右反転させるとdnegになります。

row[i],col[j],dpos[i + j],dneg[i - j + N + 1]のいずれかに襲撃フラグが立っていると、
マス目が襲撃されていることになるので、
row[i],col[j],dpos[i + j],dneg[i - j + N + 1]のすべてのフラグがオフならば、クイーンが配置できることになります。

c++による8クイーン問題のバックトラック手法の解法

では、最後にc++でのプログラムの解法を載せます。

#include <iostream>

using namespace std;

#define N 8
#define FREE -1
#define NOT_FREE 1

int row[N],col[N],dpos[2 * N - 1],dneg[2 * N - 1];

// 入力時に置かれたqueenを記録しておく配列,true = queenが置かれた
bool X[N][N];

void init(){
    for(int i = 0;i < N;++i){
        row[i] = col[i] = FREE;
    }
    for(int i = 0;i < 2 * N - 1;++i){
        dpos[i] = dneg[i] = FREE;
    }
    
    for(int i = 0;i < N;++i){
        for(int j = 0;j < N;++j){
            X[i][j] = false;
        }
    }
}

/**
 出力関数
*/
void printBoard(){
    for(int i = 0;i < N;++i){
        for(int j = 0;j < N;++j){
            // 入力時に置かれたクイーンと同じでない場合は出力しない
            if(X[i][j]){
                if(row[i] != j){
                    return;
                }
            }
        }
    }
    
    for(int i = 0;i < N;++i){
        for(int j = 0;j < N;++j){
            // rowにqueenの置かれた位置を記録している
            cout << ((row[i] == j) ?  "Q" : ".");
        }
        cout << endl;
    }
}

void recursive(int i){
    // 最後の行までクイーンを置くことができた
    if(i == N){
        printBoard();
        return;
    }
    
    for(int j = 0;j < N;++j){
        // (i,j)が他のクイーンに襲撃されている場合は飛ばす
        if(col[j]  == NOT_FREE || dpos[i + j] == NOT_FREE || dneg[i - j + N - 1] == NOT_FREE){
            continue;
        }
        // (i,j)にクイーンを配置する
        // rowにqueen置かれた位置を記録する、行番号iの何列目にqueenが置かれたのかを示す
        row[i] = j;
        col[j] = dpos[i + j] = dneg[i - j + N - 1] = NOT_FREE;
        // 次の行を試す
        recursive(i + 1);
        // 次の行がダメだった場合(i,j)に配置されたクイーンを取り除く
        row[i] = col[j] = dpos[i + j] = dneg[i - j + N - 1] = FREE;
    }
    // クイーンの配置に失敗した
}

int main(){
    init();
    
    int k;
    cin >> k;
    for(int i = 0;i < k;++i){
        int r,c;
        cin >> r >> c;
        X[r][c] = true;
    }
    
    recursive(0);
    
    return 0;
}

参考
プログラミングコンテスト攻略のためのアルゴリズムとデータ構造

2019年8月5日月曜日

プログラムの最大公約数(gcd)問題をユーグリットの互除法で解く方法

最大公約数を求める問題の解法について、解説したいと思います。
まずは、最大公約数の定義についておさらいします。

最大公約数とは

2つの整数x,yにおいて、x ÷ d,y ÷ dの余りがともに0になるdのうち最大のもの
x,yの最大公約数(great common divisor訳してgcd)と言います。

最大公約数の例

2つの整数を35,14とすると、
35の約数は(1,5,7,35)、14の約数は(1,2,7,14)
この2数の公約数(1,7)のうち最大のもの7が最大公約数になります。

最大公約数を求める例題

AIZU ONLINE JUDGEの問題を例として、最大公約数を求める問題を解いてみたいと思います。

2つの自然数 x, y を入力とし、それらの最大公約数を求めるプログラムを作成してください。

入力

x と y が1つの空白区切りで1行に与えられます。

出力

最大公約数を1行に出力してください。

単純な解法

gcd(x,y)
 n = (xとyの小さい方の数)
 for dがnから1まで
  if dがxとyの約数
   return true

この解法は、与えれれた2数x,yのうち小さな数をnとして、
dがnから1までの間で、xとyがともに約数となる数を調べるコードになります。
この単純な解法だと規定時間内にクリアすることができないので、
ユーグリットの互除法を利用して問題を解きます。

ユーグリットの互除法とは

ユーグリットの互除法とは、
$x \geq y$の時関数gcd(x,y)と関数gcd(y,xをyで割った余り)は正しい
という定理のことです。

この定理を使って最大公約数を効率的に求めることができます。

ユーグリットの互除法の使用例

35,14の2数を例として、
ユーグリットの互除法をかけてみたいと思います。

gcd(35,14)
= gcd(35,35 % 14) = gcd(35,7)
= gcd(7,35 % 7) = gcd(7,0) = 7

gcd(a,b)において、b = 0となった時のaが与えられたa,bの最大公約数になります。

ユーグリットの互除法の証明

ユーグリットの互除法が正しいことを、aとbの公約数とbとr(a % b)の公約数が
正しい
ことを利用して証明します。

まず、dをaとbの公約数とします。
dは自然数l,mを用いて、a = ld,b = mdと表せます。

a = bq + rにa = ldを代入して、
ld = bq + r -①

①にb = mdを代入して、
ld = mdq + r

rについてまとめると

r = d(l - mq)が得られます。

この式からdがrの約数だということがわかります。

また、dはa,bの公約数なので、dはbを割り切ります。
なので、dはrとbの公約数になります。

同様の方法で、d'がbとrの公約数だとすると、d'はaとbの公約数であることがわかります。

よって、aとbの公約数の集合とrとbの公約数の集合は正しく、最大公約数も同じになります。

ユーグリットの互除法の計算効率

例として、74と54にユーグリットの互除法を適応して、
計算効率をチェックしてみます。

74 = 74 % 54 + 20(r1)
54 = 54 % 20 + 14(r2)
20 = 20 % 14 + 6(r3)
14 = 14 % 6 + 2(r4)
6 = 6 % 2 + 0(r5)

ここで、ユーグリットを適応して、得られた列(r1...r5)がどのように減っていったのかを考えます。

a = bq + rとすると(0 < r < b)、$r < \frac{a}{2}$より
$r_{i + 2} = \frac{n}{2}$が成り立ちます。

これより、ユーグリットの互除法は多くとも$2log_{2}(b)$で完了するので、
O(log(b))と見積もることができます

ユーグリットの互除法の疑似アルゴリズム

ユーグリットの互除法の疑似アルゴリズムは以下のようになります。

gcd(x,y)
 if x < y
   x >= yとなるようにxとyを入れ替える

 while y > 0
  r = x % y
  x = y
  y = r

return x;

cによる最大公約数を求める問題の解法

疑似アルゴリズムのユーグリットの互除法をcに直しています。


#include <algorithm>
#include <stdio.h>
using namespace std;

int gcd(int x,int y){
    int r;
    if(x < y){
        // 入れ替え
        swap(x,y);
    }
    
    while(y > 0){
        r = x % y;
        x = y;
        y = r;
    }
    return x;
}

int main(){
    int x,y;
    scanf("%d %d",&x,&y);
    printf("%d\n",gcd(x,y));
    return 0;
}

2019年8月4日日曜日

素数と判定させるアルゴリズム

合成数(素数でない数)xは$p <= \sqrt{x}$を満たす素因子pをもつ

という性質があります。

これを利用すれば、xが素数かどうか調べるには、2からx - 1まででなく
2から$\sqrt{x}$までの値で割り切れるかどうか調べればいい
ので、
$O(\sqrt{x})$のアルゴリズムにすることができます。

このアルゴリズムを疑似言語で実装すると以下のようになります。

疑似言語による素数と判定するアルゴリズム

isPrime(x)
 if x == 2
  return true

 if x < 2 または xが偶数
  return false

 i = 3
 // 平方根以下
 while i <= xの平方根
  if xがiで割り切れる
   return false
  // 偶数はすでにチェックしているので、i + 2となる
  i = i + 2
  
 return true

2019年8月3日土曜日

最大正方形問題の動的計画法による解法と考え方

最大正方形問題の解法とその考え方について、まとめたいと思います。
AIZU ONLINE JUDGEの例題をあげます。

例題

一辺が 1 cm のタイルが、H × W 個並べられています。タイルは汚れているもの、綺麗なもののいずれかです。
綺麗なタイルのみを使ってできる正方形の面積の最大値を求めてください。

入力

H W
c1,1 c1,2 ... c1,W
c2,1 c2,2 ... c2,W
:
cH,1 cH,2 ... cH,W

1行目に2つの整数 H、W が空白区切りで与えられます。続くH 行でタイルを表す H × W 個の整数 ci,j が与えられます。
ci,j が 1 のとき汚れたタイル、0 のとき綺麗なタイルを表します。

出力

面積の最大値を1行に出力してください。

入力例

4 5
0 0 1 0 0
1 0 0 0 0
0 0 0 1 0
0 0 0 1 0

出力例

4

動的計画法を使った解法

小さな部分の問題の解を記録しておく領域をdp[H][W]として定義します。
dp[i][j]にタイル(i,j)から左上に向かってできる最大の正方形の辺長さ(タイルの数)を記録していきます。

この時の短縮方法について、図を用いて考えます。

上図のようにすでに、dpに2,1,3の正方形が入っているとします。
図で示されているように、正方形は左上に向かってできています。 これらの記録を利用して、?の部分の正方形の辺の長さを求めてみます。

dp[i][j]の値はその左上、上、左の要素の中で最も小さい値に1を加えたものになります。
図では、現在の位置(i,j)の右下とした正方形の辺の長さは、最小値である上の値1 + 1より大きくすることはできないことを示しています。

次のように、各行を左から右へ見ていく処理を、上から下へ順番に行えば、
dp[i][j]を計算する際に、左上、上、左の要素はすでに計算済みなので、これらの解を有効に利用することができます。

疑似言語による動的計画法を使った最大正方形を求めるアルゴリズム

for i = 1 to H - 1
 for j = 1 to W - 1
  if G[i][j] dirty
   dp[i][j] = 0
  else
   // 左上,上,左の最小値の+1がdp[i][j]の正方形の最大面積となる
   dp[i][j] = min(dp[i-1][j-1],min(dp[i-1][j],dp[i][j - 1])) + 1
   // 最大値を更新
   maxWidth = max(maxWidth,dp[i][j]

では、例題の解法をc++で書いてみます。

#include 
#include
using namespace std;

#define MAX 1400

int dp[MAX][MAX],G[MAX][MAX];

int getLargestSquare(int h,int w){
    int maxWidth = 0;
    for(int i = 0;i < h;++i){
        for(int j = 0;j < w;++j){
            // 1が汚れたタイルなので、面積は0に、0は綺麗なタイルなので、面積が1になる
            dp[i][j] = (G[i][j] + 1) % 2;
            // 面積が1の場合、最大面積を計算する処理がスルーされるので、面積を暫定でいれる
            maxWidth |= dp[i][j];
        }
    }

    for(int i = 1;i < h;++i){
        for(int j = 1;j < w;++j){
            // 汚れている状態
            if(G[i][j]){
                dp[i][j] = 0;
            }else{
                // 左上,上,左の最小値の+1がdp[i][j]の正方形の最大面積となる
                dp[i][j] = min(dp[i-1][j-1],min(dp[i-1][j],dp[i][j - 1])) + 1;
                maxWidth = max(maxWidth,dp[i][j]);
            }
        }
    }
    return maxWidth * maxWidth;
}

int main(){
    int h,w;
    scanf("%d %d",&h,&w);
    
    for(int i = 0;i < h;++i){
        for(int j = 0;j < w;++j){
            // 地面の状態
            scanf("%d",&G[i][j]);
        }
    }
    
    printf("%d\n",getLargestSquare(h,w));
    return 0;
}

参考
プログラミングコンテスト攻略のためのアルゴリズムとデータ構造

2019年7月30日火曜日

最長増加部分列問題(LIS)とその解法アルゴリズム

問題と解法の前に私がわからなかったいくつかの定義について解説します。

部分列とは

数列A{5,1,3,2,4}があるとすると、部分列は

{5,3,2},{1,2,4}などが相当します。

つまり、数列の先頭から部分的に数列を取り出したものになります。

なので、{1,3,5}など途中から前の数字を取り出した数列は部分列に相当しません

記号的に定義すると

数列A=a0,a1,…,an−1とすると、
0 ≤ i0 < i1 < … < ik < n かつ ai0 < ai1 < … < aik
となります。

最長増加部分列とは

先ほどのように数列Aから切り出した部分列において、数が増加し続けているものの最大のものをいいます。

数列Aを例にとると

{1,3,2}は1と3が増加しているのみなので、長さは1
{1,2,4}は1,2,4と増加しているので、長さは3となり。
{1,2,4}が最長増加部分列となります。

これらを踏まえてAIZU ONLINE JUDGEの例題を解いてみます。

最長増加部分列の例題

数列 A = {a0, a1, … , an−1} の最長増加部分列 (LIS: Longest Increasing Subsequence) の長さを求めてください。 数列 A の増加部分列は 0 ≤ i0 < i1 < … < ik < n かつ ai0 < ai1 < … < aik を満たす部分列 {ai0, ai1, … , aik} です。最長増加部分列はその中で最も k が大きいものです。

入力

1行目に数列 A の長さを示す整数 n が与えられます。続く n 行で数列の各要素 ai が与えられます。

出力

最長増加部分列の長さを1行に出力してください。

動的計画法と2分探索を用いた解法

最長増加部分列問題は、動的計画法と2分探索を用いることで、効率的に求めることができます。
以下の2つの変数を定義します。

変数名例 役割
L[n] L[i]を長さがi + 1であるような増加部分列の最後の要素の最小値とする配列
length i番目の要素までを使った最長増加部分列の長さを表す整数

処理の図解

入力をA = {4,1,6,2,8,5,6,3}に対するLの値の変化を見ていきます。

最初の施行

初めは、比較することもないので、Aの最初の入力値である4をL代入します。

A:4 L:4, , , , , , ,

2度目の施行

Aから1を選択します。
Lは増加部分列の最後の要素の最小値とする配列なので、
4と1を入れ替えます。

AをLに代入する際にAの値でLに対し2分探索を使うことで、Lの条件を満たすことができます。

A:1 L:1, , , , , , ,

3度目の施行

Aは6になります。
6は現在のLISの最後の要素よりも大きいので、
LISの最後の要素[1]に6を代入し、Lengthが1増加します。

A:6 L:1,6, , , , , ,

4度目の施行

Aは2になります。
2度目の施行と同様に2よりも大きいLの要素6と2を入れ替えます。
この時、長さは変わりません。

A:2 L:1,2, , , , , ,

最後まで施行する

以下同様の処理を繰り返します。
処理の途中結果のみ図示します。

A:8 L:1,2,8, , , , ,
A:5 L:1,2,5, , , , ,
A:7 L:1,2,5,7, , , ,
A:3 L:1,2,3,7, , , ,

最終的な、Lengthの値が最長増加部分列の長さということになります。
この動的計画法と2分探索方を用いた計算量は$O(nlog_{n})$になります。

c++での解法

最後にAIZU ONLINE JUDGEの問題を例に解法例を載せます。

#include<iostream>
#include<algorithm>

using namespace std;

#define MAX 100000

int n,A[MAX + 1],L[MAX];

/**
 最長増加部分列の長さを返す
*/
int lis(){
    // 0には初期値を
    L[0] = A[0];
    int length = 1;
    
    // 1から詮索する
    for(int i = 1;i < n;++i){
        // LISの最後の要素に追加
        if(L[length - 1] < A[i]){
            L[length++] = A[i];
        }else{
            // 2分探索を使い、選択した値をLが昇順になるように代入する
            *lower_bound(L,L + length,A[i]) = A[i];
        }
    }
    return length;
}

int main(){
    cin >> n;
    for(int i = 0;i < n;++i){
        cin >> A[i];
    }
    
    cout << lis() << endl;
    return 0;
}

参考
プログラミングコンテスト攻略のためのアルゴリズムとデータ構造

2019年7月24日水曜日

ナップザック問題の解き方

ナップザック問題とはどういう問題なのか、その解法を理解できるようにまとめたいと思います。

ナップザック問題の例題

価値がVi重さがWiであるようなN個の品物と、容量がWのナップザックがあります。
次の条件を満たすように、品物を選んでナップザックに入れます。

1.選んだ品物の価値の合計をできるだけ高くする
2.選んだ品物の重さの総和はWを超えないようにする

このような各品物を選択するかしないかの組み合わせの問題を0-1ナップザック問題といいます。

問題を考えると、N個の各品物を選択するか選択しないかの組み合わせを全て調べることになるので、
計算効率は0(N^2)ということになります。

0-1ナップザック問題は動的計画法により、0(NW)の効率で求めることができます。

0-1ナップザック問題の解法

まず、次のような役割の変数を用意します。

変数名例 役割
items[N + 1] items[i].w,items[i].vにそれぞれ重さと価値を代入するための一元配列。
構造体などを使って定義します。
C[N + 1][W + 1] i個目までの品物を考慮して、大きさWのナップザックに入れる場合の価値の合計の最大値をC[i][w]とする2次元配列

例を挙げると、C[1][3]には品物1番目におけるナップザック容量3における最大価値が入ることになり、
配列の末端であるC[N][W]には、全ての組み合わせにおける最大の価値が入ることになります。

考慮する品物の数iに対し、各iにおけるナップザックの重さwを増やしていき、
C[i][w]を更新していきます。
C[i][w]の値は

1.C[i - 1][w - 品物iの重さ] + 品物iの価値または、
2.C[i - 1][w]

の大きいほうになります。
1はこの時点で、対象の品物iを選択する。
2はこの時点で、対象の品物iを選択しないという処理になります。

また、1の品物iを選択するケースでは、品物iの重さがwを超えない場合のみ選択されます。

では、例題を元に処理に流れを見ていきます。
入力は、例題の通り

4 5
4 2
5 2
2 1
8 3

つまり、
品物の最大数N = 4、ナップザックの最大容量W = 5となります。

0-1ナップザック問題の解法

要素が0のケース

要素が0のケースを考えます。
要素が0の場合は、C[0][Wまで]の要素をすべて0に初期化します。

要素が1のケース

要素が1の場合を考えます。

品物1の重さは2なので、大きさが1のナップザックには入れることができないので、C[1][1]は0となります。 大きさが2のナップザックには、品物1を入れることができます。
ここで、選択した場合は重さが0 + 4 = 4(幅がw = 2の斜めの矢印)、
選択しない場合は0(縦の矢印)となり、大きい方の4をc[1][2]に記録します。
大きさ3,4,5のナップザックについても、同様の処理をします。

要素が2のケース

大きさが1のナップザックに、品物2は入られないので、C[2][1]は0になります。
大きさが2,3,4,5のナップザックには、品物2を入れることができます。
C[2][4]ではC[1][2](品物2の重さが2なので、4 - 2で前の品物の価値最大値をチェックしたいので、C[1][2]となります。図の斜め矢印は品物の重さを示しています。) + 品物2の価値(つまり4 + 5)または、C[1][4](価値4)つまりこの品物を選ばないケースのいずれか大きい方を選びます。
C[2][5]も同様の考え方をします。

要素が3のケース

品物3の重さは1なので、大きさが1から5のナップザックに、品物3を入れることができます。
各々の大きさで、品物を選択する場合と選択しない場合の最適なケースを考えて、選択していきます。

要素が4のケース

大きさが1,2のナップザックに、品物4をいれることはできません。
斜めの矢印は重さを意味するので、容量を超えることは、斜めの矢印が配列オーバーすることを意味します。
大きさが3,4,5のナップザックには、品物4をいれることができるので、下図のように最適な選択をします

図でまとめる

C{N][W]が最大の価値になります。
斜めの矢印が選択した品物を表しており、
C{N][W]から斜めの矢印を逆にたどっていくと、選ぶべき品物がわかるようになっています。

ということで、品物の選択状況を別の2次元配列[i][w]に記憶しておけば、
最適解における品物の組み合わせを復元することができます。
2次元配列をG[i][w]とすると、品物iが選択するときをDIAGONAL、選択されなかった時をTOPとして、
変数定義をすることが一般的なようです(boolでももちろん可)。

例題をc++を使って解く

例題をc++を使って、解説を入れつつ解いてみたいと思います。

まず、必要な変数を定義します。

#define NMAX 105
#define WMAX 10005
// 選択する場合
#define DIAGONAL 1
// 選択しない場合
#define TOP 0

// 価値と重さを持った構造体を定義
struct Item {
    int value,weight;
};
// 品物の数(入力により決められる)
int N;
// ナップザックの容量(入力により決められる)
int W;
// 入力された品物の価値と重さを保存する配列
Item items[NMAX + 1];
// 各々に価値の最適解を入れておく配列
int C[NMAX + 1][WMAX + 1];
// 各々の品物を選択するか選択しないかを入れておく配列
int G[NMAX + 1][WMAX + 1];

続いて、図解で示した処理を行う関数を定義します。 コメントで処理内容の詳細を書きました

/**
    0-1ナップザック問題の解法を行う
    最大値を返す。
*/
int compute(){
    // 配列要素0を0で初期化(図参照)
    for(int w = 0;w <= W;++w){
        C[0][w] = 0;
        G[0][w] = DIAGONAL;
    }
    
    // C[1 to N)を0で初期化(図参照)
    for(int i = 1;i <= N;++i) C[i][0] = 0;
    
    // 図の要素1から要素の末尾までの処理を行う
    for(int i = 1;i <= N;++i){
        for(int w = 1;w <= W;++w){
            // 選択しなかったケースを代入しておく
            // 選択しなかったケースは図のように[i - 1]を重さを入れておく
            C[i][w] = C[i - 1][w];
            // Gは選択するかしないかを格納する配列、格納しないケースであるTOPを入れる
            G[i][w] = TOP;
            // 現時点でのナップザックの重さが対象の品物の重さ以上でないと、その品物を選択することができない
            if(items[i].weight > w) continue;
            // 左辺は対象の品物の価値と品物の重さを引いたC配列の最適解を加算したものを示しています。
            // 右辺はC[i - 1][w]つまり選択しなかったケースです。
            // この大小をチェック、つまりは最適解のチェックをしています。
            if(items[i].value + C[i - 1][w - items[i].weight] > C[i - 1][w]){
                // 更新する
                C[i][w] = items[i].value + C[i - 1][w - items[i].weight];
                G[i][w] = DIAGONAL;
            }
        }
    }
    
    // C[N][W]が最大となる
    return C[N][W];
}

問題には関係ありませんが、選んだ品物を順番に表示するには、 以下のように処理を書きます。

/**
 選んだ品物を順番に表示する
*/
void printSelection(){
    vector selection;
    for(int i = N, w = W;i >= 1;--i){
        if(G[i][w] == DIAGONAL){
            selection.push_back(i);
            // 選んだ品物の重量を引く、斜め矢印の処理にあたる
            w -= items[i].weight;
        }
    }
    
    reverse(selection.begin(),selection.end());
    
    for (int i = 0;i < selection.size();++i) {
        printf("%d ",selection.at(i));
    }
}

プログラムの全体は以下のようになります。

#include<iostream>
#include<vector>
#include<algorithm>

using namespace std;

#define NMAX 105
#define WMAX 10005
// 選択する場合
#define DIAGONAL 1
// 選択しない場合
#define TOP 0

// 価値と重さを持った構造体を定義
struct Item {
    int value,weight;
};
// 品物の数(入力により決められる)
int N;
// ナップザックの容量(入力により決められる)
int W;
// 入力された品物の価値と重さを保存する配列
Item items[NMAX + 1];
// 各々に価値の最適解を入れておく配列
int C[NMAX + 1][WMAX + 1];
// 各々の品物を選択するか選択しないかを入れておく配列
int G[NMAX + 1][WMAX + 1];

/**
    0-1ナップザック問題の解法を行う
    最大値を返す。
*/
int compute(){
    // 配列要素0を0で初期化(図参照)
    for(int w = 0;w <= W;++w){
        C[0][w] = 0;
        G[0][w] = DIAGONAL;
    }
    
    // C[1 to N)を0で初期化(図参照)
    for(int i = 1;i <= N;++i) C[i][0] = 0;
    
    // 図の要素1から要素の末尾までの処理を行う
    for(int i = 1;i <= N;++i){
        for(int w = 1;w <= W;++w){
            // 選択しなかったケースを代入しておく
            // 選択しなかったケースは図のように[i - 1]を重さを入れておく
            C[i][w] = C[i - 1][w];
            // Gは選択するかしないかを格納する配列、格納しないケースであるTOPを入れる
            G[i][w] = TOP;
            // 現時点でのナップザックの重さが対象の品物の重さ以上でないと、その品物を選択することができない
            if(items[i].weight > w) continue;
            // 左辺は対象の品物の価値と品物の重さを引いたC配列の最適解を加算したものを示しています。
            // 右辺はC[i - 1][w]つまり選択しなかったケースです。
            // この大小をチェック、つまりは最適解のチェックをしています。
            if(items[i].value + C[i - 1][w - items[i].weight] > C[i - 1][w]){
                // 更新する
                C[i][w] = items[i].value + C[i - 1][w - items[i].weight];
                G[i][w] = DIAGONAL;
            }
        }
    }
    
    // C[N][W]が最大となる
    return C[N][W];
}

int main(){
    int maxValue;
    // 入力処理
    cin >> N >> W;
    for(int i = 1;i <= N;++i){
        cin >> items[i].value >> items[i].weight;
    }
    maxValue = compute();
    
    cout << maxValue << endl;
    return 0;
}

参考
プログラミングコンテスト攻略のためのアルゴリズムとデータ構造