裏 RjpWiki

文字通り,RjpWiki の裏を行きます
R プログラム コンピュータ・サイエンス 統計学

ソートされないカード

2016年12月06日 | ブログラミング

ソートされないカード

締め切りが 2016/12/06 10:00 AM なので,その 1 分後に投稿されるように予約

設問
1〜n までの整数が1つずつ書かれている n 枚のカードを横一列に並べます。
カードを左から順に一枚ずつ見て、書かれている数字が i なら、左から i 番目のカードと交換する、
という作業を右端のカードまで繰り返します。

例えば、3, 2, 5, 4, 1の順に並んでいると、
最初のカードは 3 なので3番目のカードである 5 と交換し、5, 2, 3, 4, 1となります。
次のカードは 2 なので2番目のカードと交換(交換は発生しない)、
その次のカードは 3 なので3番目のカードと交換(交換は発生しない)、
その次のカードは 4 なので4番目のカードと交換(交換は発生しない)、
その次のカードは 1 なので1番目のカードである 5 と交換し、
1, 2, 3, 4, 5となり、昇順に並べ替えられます。
※カードは常に左から順番に見ていきます

右端まで処理した結果、書かれている数字が昇順に並ばない初期配置が何通りあるかを求めます。

標準入力から n が与えられるとき、書かれている数字が昇順に並ばない初期配置が何通りあるかを求め、
標準出力に出力してください。
例えば、n = 4 のとき、24通り中、以下の3通りがありますので、サンプルのように出力します。

2, 3, 4, 1
3, 4, 2, 1
4, 3, 1, 2
(n は最大で8とします。余裕がある人は、もう少し大きな n についても考えてみてください。)

【入出力サンプル】
標準入力
4

標準出力
3

============================

ブルートフォースでは,n=8 の場合に実行時間制限(1秒)に引っかかるのだが,細々したところをチューンアップして,なんとかねじ込んだ。

f = function(n) {
    z = matrix(1)
    if (n == 1) {
        cat(0)
    } else {
        for (i in 2:n) {
            x = cbind(z, i)
            a = c(1:i, 1:(i - 1))
            z = matrix(0, ncol = ncol(x), nrow = i * nrow(x))
            z[1:nrow(x), ] = x
            for (j in 2:i - 1) {
                z[j * nrow(x) + 1:nrow(x), ] = x[, a[1:i + j]]
            }
        }
        cat(sum(apply(z, 1, function(x) {
            for (i in seq_along(x)) {
                if (i != x[i]) {
                    t = x[i]
                    x[i] = x[t]
                    x[t] = t
                }
            }
            any(x != 1:n)
        })))
    }
}

f(3) # 0
f(4) # 3
f(5) # 35
f(6) # 325
f(7) # 2989
f(8) # 28630

コメント
この記事をはてなブックマークに追加

「ディビジョン・サム」問題

2016年12月01日 | ブログラミング

「ディビジョン・サム」問題

締め切りが 2016/12/01 10:00 AM なので,その 1 分後に投稿されるように予約

設問
自然数 n に対し、(n!)^n (つまり n の階乗の n 乗)の約数の和を F(n) と定義します。
例えば、F(2) = 7 です。(2!)^2 = 4 で、4 の約数は 1, 2, 4 だからです。
また、F(3) = 600 です。(3!)^3 = 216 で、216 には 16 個の約数があり、それらの和は 600 です。
同様に、F(4) = 991111となります。
標準入力から、自然数 n(1 ≦ n ≦ 103)が与えられるとき、
標準出力に F(n) を 1000003 で割った余り を出力するプログラムを書いてください。
たとえば、
F(12) を1000003で割った余りは845575となります。

===============================

R で書くと F(997) が時間オーバーになった
AWK で書き,bash で gawk を起動すれば,制限時間内で答えが出た

M  = 1000003
divisor = function(n) {
    if (n%%2 == 0)
        return(2)
    else if (n%%3 == 0)
        return(3)
    maxitr = as.integer(sqrt(n))
    i = 1
    repeat {
        i = i + 4
        if (i > maxitr)
            return(n)
        else if (n%%i == 0)
            return(i)
        i = i + 2
        if (i > maxitr)
            return(n)
        else if (n%%i == 0)
            return(i)
    }
}

factorization = function(n) {
    result = NULL
    repeat {
        div = divisor(n)
        result = c(result, div)
        n = n/div
        if (n == 1)
            return(result)
    }
}

prod = function(b, a) {
    result = term = 1
    for (i in 1:a) {
        term = (term*b) %% M
        result = (result+term)
    }
    result %% M
}

F = function(n) {
    result = NULL
    for (i in 2:n) {
        result = c(result, factorization(i))
    }
    a = n*table(result)
    b = as.integer(names(a))
    ans = 1
    for (i in seq_along(a)) {
        ans = (ans*prod(b[i], a[i])) %% M
    }
    cat(ans %% M)
}
# F(4) # 991111
# F(12) # 845575
# F(11) # 687005
# F(99) # 820272
# F(264) # 398251
# F(610) # 823587
# F(912) # 210178
# F(997) # 772380

==========================

awk

function divisor(n,     i, maxitr) {
  if (n%2 == 0) return 2
  else if (n%3 == 0) return 3
  maxitr = int(sqrt(n))
  i = 1
  for (;;) {
    i = i + 4
    if (i > maxitr) return n
    else if (n%i == 0) return i
    i = i + 2
    if (i > maxitr) return n
    else if (n%i == 0) return i
  }
}

function factorization(n, result,    div) {
  for (;;) {
    div = divisor(n)
    result[div]++
    n = n/div
    if (n == 1) return
  }
}

function prod(b, a,     result, term, i, M) {
  M = 1000003
  result = term = 1
  for (i = 1; a >= i; i++) {
    term = term*b % M
    result = result+term
  }
 return result % M
}

function F(n,       i, res, fac, result, ans, b, M) {
  M = 1000003
  for (i = 2; n >= i; i++) {
    factorization(i, res)
    for (fac in res) {
      result[fac] += res[fac]
      delete res[fac]
    }
  }
  ans = 1
  for (b in result) {
    ans = (ans*prod(b, n*result[b])) % M
  }
  print ans % M
}

コメント
この記事をはてなブックマークに追加

不要なスペース除去

2016年12月01日 | ブログラミング

締め切りが 2016/12/01 10:00AM なので,その 1 分後に投稿するように予約

求められるプログラムの前提条件は、以下の通りとなります。

    標準入力から、英数字、スペース、句読点(後述)のいずれかのみで構成された英文が送られる
    文字列の長さは80以下とする
    文字列の先頭および末尾から続くスペースは除去の対象する
    先頭および末尾を処理した文字列の中で、2文字以上連続するスペースは、1スペースを残し、残るスペースは除去の対象とする
    句読点はピリオド(.)、カンマ(,)、疑問符(?)のいずれかとし、句読点の直前にあるスペースは除去対象とする
    上記ルールにもとづき、スペースを除去した文字列を標準出力に返すこと

以下、整形例となります。

【入出力サンプル】
標準入力
  Hello  ,  CodeIQ.  

標準出力
Hello, CodeIQ.

==========

「○○となります」ってのは,気味悪いから,止しにしてくれ。
なお,テスト入力例が,プログラムの仕様をチェックするために十分なものでないのが納得いかねえ。

f = function(s) {
    s = sub("^ +", "", s)
    s = sub(" +$", "", s)
    s = gsub("  +", " ", s)
    s = gsub(" +([,.?])", "\\1", s)
    cat(s)
}

コメント
この記事をはてなブックマークに追加

ロールシャッハ

2016年11月25日 | ブログラミング

図版を公開するのは止めろとかいう無駄な圧力があったり,

解釈もでたらめだという話が「心理テスト」はウソでした。 受けたみんなが馬鹿を見た」 単行本 – 2005/3/30  村上 宣寛 (著)  にも書いてあった。

 

コメント
この記事をはてなブックマークに追加

Base64で反転

2016年11月22日 | ブログラミング

Base64で反転

締め切りが 2016/11/22 10:00 AM なので,その 1 分後に投稿されるように予約

設問
A-Zとa-zの52文字から構成される、長さが 3n の文字列があります。
これをASCIIコードからBase64にエンコードし、左右反転します。
さらにBase64からデコードしたときに、元の文字列と同じになるもののうち、
元の文字列に含まれる文字が n 種類のものがいくつあるかを出力してください。

例えば、n = 1 のとき、「TQU」という文字列はエンコードすると「VFFV」となり、
左右反転してデコードすると「TQU」に戻ります。
ただ、この場合は「T」「Q」「U」という3種類の文字を使用しています。

同様に、「DQQ」「fYY」は2種類の文字を使用しています。
n = 1 のときは「UUU」の1つだけですので、1を出力します。
なお、n は5以下の整数とします。

【入出力サンプル】
標準入力
1

標準出力
1

===================

R で書いたら, n = 5 のときに制限時間オーバー

f = function(N) {
  tbl = c(LETTERS, letters)
  int.tbl = sapply(tbl, utf8ToInt)
  tbl1 = c("a", "A", "b", "B", "d", "D", "e", "E", "f", "F", "h", "H", "i", "I", "j", "J", "L", "M",
    "N", "P", "Q", "R", "T", "U", "V", "X", "Y", "Z")
  tbl2 = c("Q", "U", "Y")
  tbl3 = c("P", "Q", "R", "S", "T", "U", "V", "X", "Y", "Z")
  a0 = as.matrix(expand.grid(tbl1, tbl2, tbl3))
  a = as.matrix(expand.grid(sapply(tbl1, utf8ToInt), sapply(tbl2, utf8ToInt), sapply(tbl3, utf8ToInt)))
  intToBin = matrix(0L, 8, max(int.tbl))
  for (i in int.tbl) {
    intToBin[, i] = as.integer(intToBits(i))[8:1]
  }
  intToBin2 = matrix(0L, 6, 52)

  for (i in 1:52) {
    intToBin2[, i] = as.integer(intToBits(i - 1))[6:1]
  }

  # アルファベット3文字の組の2進表現(3バイト;24ビット) 24×140608 (=52^3)
  b = rbind(intToBin[, a[, 1]], intToBin[, a[, 2]], intToBin[, a[, 3]])

  # エンコーディングされたアルファベット4文字(6ビット)の組として読む
  w = 2^(5:0)
  z1 = w %*% b[1:6, ] # 0 ~ 51 ならアルファベットに変換される
  z2 = w %*% b[7:12, ]
  z3 = w %*% b[13:18, ]
  z4 = w %*% b[19:24, ]
  valid = (52 > z2) & (52 > z3) & (52 > z4) # z1 はチェック無用
  z1 = z1[valid]
  z2 = z2[valid]
  z3 = z3[valid]
  z4 = z4[valid]
  a = a0[valid, ] # アルファベット3文字の組 82800×3
  w1 = intToBin2[, z1 + 1] # 4文字の組からデコードしてアルファベット3文字の組を作る
  w2 = intToBin2[, z2 + 1]
  w3 = intToBin2[, z3 + 1]
  w4 = intToBin2[, z4 + 1]
  y = rbind(w4, w3, w2, w1)
  w8 = 2^(0:7)
  range = int.tbl # sapply(tbl, utf8ToInt)
  y1 = w8 %*% y[8:1, ]
  y2 = w8 %*% y[16:9, ]
  y3 = w8 %*% y[24:17, ]
  valid2 = (y1 %in% range) & (y2 %in% range) & (y3 %in% range)
  s1 = sapply(y1[valid2], intToUtf8)
  s2 = sapply(y2[valid2], intToUtf8)
  s3 = sapply(y3[valid2], intToUtf8)
  final = cbind(a[valid2, ], s1, s2, s3) # 右3列は,左3列と逆順のエンコード文字列を生成する
  # たとえば,"AQQ" は "QVFR" となり,"DUP" は "RFVQ" となる
 
  flag = rowSums(final[, 1:3] == final[, 4:6]) == 3
  sym = final[flag, 1:3]
  rev = final[!flag, ]

  n = nrow(final)
  n.sym = nrow(sym)
  result = 0

  unique2 = function(x) length(.Internal(unique(x, FALSE, FALSE, NA)))

  if (N == 1) { # 1 ==> 1
    f = apply(sym, 1, function(x) length(table(x)) == 1)
    result = sum(f)
  } else if (N == 2) { # 2 ==> 2
    f = apply(final, 1, function(x) length(table(x)) == 2)
    result = sum(f)
  } else if (N == 3) { # 3 ==> 127
    for (i in 1:n) {
      if (unique2(final[i, ]) > 3)
        next
      for (j in 1:n.sym) {
        result = result + (unique2(c(final[i, ], sym[j, ])) == 3)
      }
    }
  } else if (N == 4) { # 4 ==> 1536
    for (i1 in 1:n) {
      if (unique2(final[i1, ]) > 4)
        next
      for (i2 in i1:n) {
        result = result + (unique2(c(final[i1, ], final[i2, ])) == 4) * ifelse(i1 == i2, 1, 2)
      }
    }
  } else if (N == 5) { # 5 ==> 52500
    f = 5 >= apply(final, 1, function(x) length(unique(x)))
    final = final[f, ]
    n = nrow(final)
    for (i1 in 1:n) {
      if (unique2(final[i1, ]) > 5)
        next
      for (i2 in i1:n) {
        part2 = c(final[i1, ], final[i2, ])
        if (unique2(part2) > 5)
          next
        for (j in 1:n.sym) {
          result = result + (unique2(c(part2, sym[j, ])) == 5) * ifelse(i1 == i2, 1, 2)
        }
      }
    }
  }
  cat(result)
}

f(1) # 1
f(2) # 2
f(3) # 127
f(4) # 1536
f(5) # 52500

コメント
この記事をはてなブックマークに追加

三角形と点の関係

2016年11月18日 | ブログラミング

締め切りが 11/18 10:00AM なので,その 1 分後に投稿するように予約する

想定時間「10分」だと...いつもながら,大幅に超過したわ。

三角形と点の関係

【概要】

三角形と点があります。
三角形と点の位置関係は、下表のいずれかになります:

記号     説明
A     点は、三角形の内側にあります。
B     点は、三角形の辺上にありますが、頂点上ではありません。
C     点は、三角形の頂点にあります。
D     点は、三角形の外側にあります。

どの関係になるのかを求めるプログラムを書いて下さい。

【入出力】

入力は
(1,1)(10,1)(1,10)/(2,2) (1,2)(1,1)(2,1)/(12,34) (1,3)(3,1)(4,4)/(2,2)
こんな感じです。

三角形と点の組が空白区切りで並んでいます。
スラッシュの前が三角形です。() 内に、頂点のx座標y座標がコンマ区切りで並んでいます。スラッシュの後は点の座標です。

出力は、各組が前述の A〜D のどれに該当するのかを区切り文字なしで並べたものです。

先程の例だと

    最初の組(1,1)(10,1)(1,10)/(2,2)の関係は、A。
    二番目の組(1,2)(1,1)(2,1)/(12,34)の関係は、D。
    最後の組(1,3)(3,1)(4,4)/(2,2)の関係は、B。

ということで、ADBを出力すれば正解となります。
末尾の改行はあってもなくても構いません。

【例】
入力     出力
(1,1)(10,1)(1,10)/(2,2) (1,2)(1,1)(2,1)/(12,34) (1,3)(3,1)(4,4)/(2,2)     ADB
(10,200)(10,10)(123,10)/(4,5) (12,35)(31,12)(41,44)/(41,44)     DC

【補足】

    座標は、いずれも正の整数で、百万を超えることはありません。
    ひとつの入力に含まれる 三角形-点 の組は 10以下です。
    いずれの三角形も、ゼロより大きな面積を持ちます。
    不正な入力に対処する必要はありません。


点 (a, b, c) の面積と,(a, b, z), (b, c, z), (c, a, z) の面積を求め,後三者の面積の和と前者の面積の関係でやろうかなと思ったけど,テストデータが意地悪(座標値がでかい)なので,別の,プリミティブな方法で書いた(かったるい)

slope = function(x1, y1, x2, y2) (y2 - y1)/(x2 - x1)

decomp = function(s) { # 入力行の分解
    t = unlist(strsplit(s, " "))
    sapply(t, function(u) {
        w = gsub("\\(", "", gsub(")", ",", sub("/", "", u)))
        as.numeric(unlist(strsplit(w, ",")))
    })
}

disc = function(p) { # 位置関係の判定
    x = p[1:3 * 2 - 1]
    y = p[1:3 * 2]
    o = order(x)
    x = c(x[o], p[7])
    y = c(y[o], p[8])
    a.x = x[1]
    a.y = y[1]
    b.x = x[2]
    b.y = y[2]
    c.x = x[3]
    c.y = y[3]
    z.x = x[4]
    z.y = y[4]
    eq.x = c(a.x, b.x, c.x) %in% z.x
    eq.y = c(a.y, b.y, c.y) %in% z.y
    if (any(eq.x & eq.y))
        ans = "C"
    else if (sum(eq.x) == 2 || sum(eq.y) == 2) {
        ans = "B"
    } else if (z.x = max(x) || z.y = max(y)) {
        ans = "D"
    } else if (z.x >= a.x && z.x
        slope.az = slope(a.x, a.y, z.x, z.y)
        slope.ab = slope(a.x, a.y, b.x, b.y)
        slope.ac = slope(a.x, a.y, c.x, c.y)
        if (slope.az == slope.ab || slope.az == slope.ac) {
            ans = "B"
        } else if (slope.az > max(slope.ab, slope.ac) || slope.az < min(slope.ab, slope.ac)) {
            ans = "D"
        } else ans = "A"
    } else if (z.x > b.x && z.x
        slope.cz = slope(c.x, c.y, z.x, z.y)
        slope.ca = slope(c.x, c.y, a.x, a.y)
        slope.cb = slope(c.x, c.y, b.x, b.y)
        if (slope.cz == slope.ca || slope.cz == slope.cb) {
            ans = "B"
        } else if (slope.cz > max(slope.ca, slope.cb) || slope.cz < min(slope.ca, slope.cb)) {
            ans = "D"
        } else ans = "A"
    } else ans = "x"
    ans
}

f = function(s) {
    p = decomp(s) # 入力行の分解 >> 4 点の x, y 座標を要素数 8 の列ベクトルにする行列
    cat(paste(apply(p, 2, disc), collapse = "")) # 位置関係の判定
}

f("(1,1)(10,1)(1,10)/(2,2) (1,2)(1,1)(2,1)/(12,34) (1,3)(3,1)(4,4)/(2,2)") # ADB

f("(10,200)(10,10)(123,10)/(4,5) (12,35)(31,12)(41,44)/(41,44)") # DC

f("(16,13)(1,11)(13,18)/(12,16) (2,16)(2,5)(17,10)/(8,7) (11,18)(4,4)(16,15)/(6,6) (5,9)(1,15)(6,13)/(9,18) (12,18)(7,12)(16,10)/(14,9) (6,18)(12,5)(8,12)/(6,18) (18,11)(18,6)(5,6)/(18,10) (8,2)(19,14)(5,11)/(8,2)") # ABADDCBC

f("(6,2)(26,15)(39,25)/(9,37) (46,39)(29,40)(49,1)/(39,27) (7,15)(13,47)(30,36)/(10,31) (1,31)(6,43)(25,37)/(25,37) (16,43)(31,16)(35,44)/(23,34) (49,36)(26,14)(4,1)/(9,15) (7,45)(41,13)(31,33)/(33,29)") # DABCADB

f("(63944,492389)(380462,882527)(496266,305904)/(63944,492389) (399396,945126)(3461,410781)(770003,22887)/(161835,624519) (11251,380400)(389242,357769)(685090,412104)/(460477,401536) (810430,271880)(924670,512828)(308863,176996)/(810430,271880)") # CBBC

f("(622671,535852)(364079,448665)(10452,792969)/(364079,448665) (41,20)(57,9)(79,50)/(52,28) (368208,870739)(996210,589804)(649917,745441)/(577542,777094) (38,60)(18,16)(31,78)/(27,40) (77,20)(54,83)(9,12)/(42,59) (279707,902871)(849496,790389)(763927,751666)/(570239,812148) (946989,831217)(602996,336137)(34367,168235)/(34367,168235)") # CABAABC

f("(178948,265045)(931263,891839)(781020,277145)/(329466,268070) (338,475)(663,835)(259,960)/(317,642) (289859,521270)(664758,949028)(854363,96377)/(478027,379639) (960,678)(185,953)(968,166)/(766,687) (481,661)(941,606)(424,399)/(505,539) (894179,443799)(722784,445061)(290417,278627)/(894179,443799) (796185,261874)(373280,282588)(13608,764317)/(13608,764317)") # BABAACC

f("(433423,181160)(832013,121400)(359517,786345)/(473282,175184) (485,8690)(4722,1849)(7785,8488)/(6430,7329) (818475,458807)(839588,892108)(379565,202673)/(818475,458807) (131839,712754)(48808,955181)(506046,158685)/(506046,158685) (8435,3330)(7171,7896)(8676,9365)/(8146,7191) (9937,4933)(3373,1020)(4835,9395)/(6082,7342) (358973,847370)(876921,331437)(247831,429607)/(436558,400156)") # BACCAAB

f("(62493,977407)(518706,994645)(299775,962634)/(366635,988899) (94459,94078)(97847,48616)(46433,75423)/(94597,79113) (828610,187056)(832354,79565)(95750,759545)/(279901,589550) (95604,54883)(67680,26554)(45312,90675)/(69873,36818) (807133,92093)(71884,217935)(554759,810479)/(71884,217935) (83067,73903)(92600,11864)(23347,94602)/(65277,73576) (140792,148212)(203528,414277)(11301,473993)/(140792,148212)") # BABACAC

f("(158839,197999)(608294,629979)(117450,803952)/(132343,724233) (16656,583111)(421768,921769)(751443,911561)/(16656,583111) (913703,284010)(698577,8906)(349662,568541)/(913703,284010) (355069,187114)(879282,384550)(537733,854548)/(575674,616336) (34536,772021)(760160,148154)(361255,817079)/(441036,683294) (439610,950121)(646821,369587)(138285,752935)/(390271,860139) (393398,255982)(711129,576840)(355440,876162)/(671608,610098)") # ACCABAB

f("(159635,876685)(888617,781119)(480723,957682)/(480723,957682) (93588,685668)(64882,993088)(162406,77103)/(106678,600523) (780162,833747)(757864,76945)(435151,731142)/(593723,736001) (434282,643984)(445591,676783)(753853,60119)/(753853,60119) (275601,291635)(778880,26144)(97423,935493)/(584178,285958) (961646,450205)(291968,836498)(87569,238183)/(207861,466806) (820255,200270)(644706,594803)(520996,202782)/(626434,319571)") # CBACBAA

f("(688483,858004)(835886,244022)(553322,192384)/(835886,244022) (51419,482928)(514395,53297)(344931,940813)/(387297,718934) (721856,686900)(856884,522826)(314664,338075)/(392124,364468) (579199,994594)(97241,936422)(502551,112180)/(330592,596296) (812841,538190)(845697,938354)(122714,149293)/(122714,149293) (105757,293420)(534586,548742)(597299,773659)/(501583,662648) (599633,724163)(694400,337219)(92480,632980)/(548438,710269)") # CBBACAA

f("(892402,565227)(480608,382305)(818167,303522)/(783524,181393) (972526,151377)(225533,281022)(270554,240711)/(255547,254148) (948978,580225)(22448,737666)(299192,435497)/(391440,334774) (623398,895656)(827192,867194)(727880,126074)/(723742,95194) (761473,966845)(349291,914983)(193691,32011)/(271491,473497) (581768,864305)(802761,210568)(224151,90772)/(706326,190602) (481556,71688)(511498,710360)(181945,471407)/(72094,391756) (210710,232479)(337274,561336)(449777,309744)/(606788,360489) (396169,793629)(680416,666213)(418241,127653)/(575546,450789) (959550,369216)(998967,880416)(154188,493257)/(422642,451910)") # DBDDBBDDBB

f("(461484,305711)(435830,117953)(582877,502129)/(536509,427104) (21441,229521)(73349,117583)(587578,949623)/(391160,631812) (351400,100879)(426425,222272)(203246,620105)/(397768,175904) (192849,99802)(707078,931842)(985751,877379)/(510660,614031) (485816,377542)(476282,467118)(794093,981347)/(672700,784929) (494245,336449)(297827,18638)(544752,124606)/(419220,215056) (162698,379593)(237723,500986)(827620,763861)/(209066,454618) (634932,955311)(120703,123271)(351149,515585)/(438514,637500) (149490,780650)(70083,25940)(116451,100965)/(98740,72308) (459254,720823)(377165,233374)(348508,187006)/(366219,215663)") # AAADADDADA

コメント
この記事をはてなブックマークに追加

「デジタル・ルート」問題

2016年11月03日 | ブログラミング

「デジタル・ルート」問題

締め切りが 2016/11/03 10:00 AM なので,その 1 分後に投稿されるように予約

ある自然数に対し、各桁の数字をすべて足し合わせるという変換を行います。
例えば、199 にこの変換を行うと、1+9+9=19 という数になります。

さらに、得られた数に同様の変換を行います。最終的に 1 桁の数になるまで変換を繰り返し行います。
199 だと、1 回目の変換で 19 になり、2 回目の変換で 10 になり、3 回目の変換で 1 になります。

自然数 k に対し、上の変換を繰り返したときに最終的に 1 桁の数になるまでに必要な変換の回数を F(k) と定義します。
例えば F(199) = 3 です。同様に、F(9) = 0, F(26) = 1, F(888888) = 3 です。

自然数 n に対し、1 ≦ k ≦ n となる全ての k に対する F(k) の和を G(n) と定義します。
例えば G(9) = 0,G(20) = 12,G(149) = 200, G(9876) = 19951 となることが確かめられます。

標準入力から、自然数 n(1 ≦ n ≦ 108)が与えられます。
標準出力に G(n) の値を出力するプログラムを書いてください。

==========================================

単純なブルートフォースによるプログラム

tbl = integer(1e8)
F = function(n) {
    in.n = n
    count = 0
    repeat {
        if (10 > n) {
            return(count)
        }
        count = count+1
        n = sum(as.integer(unlist(strsplit(as.character(n), ""))))
        if (tbl[n]) {
            count = count+tbl[n]
            count ->> tbl[in.n]
            return(count)
        }
    }
}
# F(199) # 3
# F(9) # 0
# F(26) # 1
# F(888888) # 3
G = function(k) {
    sum(sapply(1:k, F))
}
#G(48) # 48
#G(100) # 136
#G(5432) # 10539
#G(54321) # 113023
#G(90817263)
#G(99003157)

=====================================

しかし,そんなプログラムでは歯の立たないことが多いので,工夫する

seq.gen2 = function(m) {
    z = vector("list", m)
    z[[1]] = x = rep(1, 10)
    for (i in 2:m) {
        y = c(x, rep(0, 9))
        for (j in 1:9) {
            y[j+seq_along(x)] = y[j+seq_along(x)]+x
        }
        z[[i]] = y
        x = y
    }
    z
}
seq = seq.gen2(8)

リストである seq の要素は,
[[1]] 1 1 1 1 1 1 1 1 1 1
[[2]] 1 2 3 4 5 6 7 8 9 10 9 8 7 6 5 4 3 2 1
[[3]] 1 3 6 10 15 21 28 36 45 55 63 69 73 75 75 73 69 63 55 45 36 28 21 15 10 6 3 1
[[4]] 1 4 10 20 35 56 84 120 165 220 282 348 415 480 540 592 633 660 670 660 633 592 540 480 415 348 282 220 165 120 84 56 35 20 10 4 1
[[5]] 1 5 15 35 70 126 210 330 495 715 996 1340 1745 2205 2710 3246 3795 4335 4840 5280 5631 5875 6000 6000 5875 5631 5280 4840 4335 3795 3246 2710 2205 1745 1340 996 715 495 330 210 126 70 35 15 5 1
  :
[[8]] 1 8 36 120 330 792 1716 3432 ...

G(342) は,

======================================================================
*a    0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20
======================================================================
      1  2  3  4  5  6  7  8  9 10  9  8  7  6  5  4  3  2  1        
         1  2  3  4  5  6  7  8  9 10  9  8  7  6  5  4  3  2  1    
            1  2  3  4  5  6  7  8  9 10  9  8  7  6  5  4  3  2  1
               1  1  1  1  1  1  1  1  1  1                                
                  1  1  1  1  1  1  1  1  1  1                            
                     1  1  1  1  1  1  1  1  1  1                        
                        1  1  1  1  1  1  1  1  1  1                    
+                          1  1  1                                            
======================================================================
*b    1  3  6 10 14 18 22 26 29 32 32 31 28 24 20 16 12  9  6  3  1

======================================================================
を用いて計算する。
*a は,G の引数となる数の各桁を足しあわせたもの。たとえば,342 なら 3+4+2 = 9。この結果は G(0) 〜 G(99999999) に対しては,0 ~ 72(= 8*9) に落ち着く。
たとえば,この数値が 67 である場合は,さらに 6+7 = 13, 1+3 = 4 となるので,実際は全体で 3 回の桁の足し算がある。
*b は *a が何回出現するかの頻度(カウント数)。
求める答えは sum((y[num]+1)*tbl)-10 である。

g = function(n) {
    y = rep(c(0,1,2,1,2,1,2,1,2,1,2,1,2,1), c(10,9,1,8,2,7,3,6,4,5,5,4,6,3))
    m = as.integer(unlist(strsplit(as.character(n), "")))
    tbl = integer(100)
    begin = 0
    for (i in (length(m)-1):1) {
        x = seq[[i]]
        for (j in seq_len(m[length(m)-i])) {
            tbl[begin+seq_along(x)] = tbl[begin+seq_along(x)]+x
            begin = begin+1
        }
    }
    x = rep(1, m[length(m)]+1)
    tbl[begin+seq_along(x)] = tbl[begin+seq_along(x)]+x
    tbl = tbl[tbl != 0]
    num = seq_along(tbl)
    cat(sum((y[num]+1)*tbl)-10)
}

G(48) # 48
g(48)

G(100) # 136
g(100)

G(5432) # 10539
g(5432)

G(54321) # 113023
g(54321)

g(90817263) # 207841726

g(99003157) # 227032799

コメント
この記事をはてなブックマークに追加

金山発掘(馬鹿馬鹿し)

2016年11月01日 | ブログラミング

締め切りが 11/01 12:00 AM なので,その 1 分後に登録されるように予約する

問題

あなたは誕生日プレゼントとして、土地をもらえることになりました。
もらえる土地は、5×5マスの正方形の土地です。
世界は10×10のマスでできており、いくつかのマスには金山があります。金山の数が最大になる5×5マスは、1つしかありません。
最も多くの金山をゲットできるように、もらえる土地を探索するプログラムを書いてください。



以下、入力の例です。左上を {"x":0,"y":0} として、右下を {"x":9,"y":9} とします。「x」はx座標、「y」はy座標です。通常の土地は「w」の文字、金山のある土地は「G」の文字です。
これらは標準入力から、改行で区切られた文字として渡されます。

wwGwwwwwGG
Gwwwwwwwww
wwwwwwwwww
Gwwwwwwwww
wwwwGwwGww
wGwwwwwwww
wwwGGwwwww
wwwwwwGwww
wwwwGGwwww
GwwwGGwGwG

以下、出力の例です。最も多くの金山が得られる土地の左上の座標を、「{"x":3,"y":5,"g":8}」のように答えます。「x」はx座標、「y」はy座標、「g」は5×5マスの土地に含まれる金山の数です。
答えは、以下のように標準出力に出力してください。
{"x":3,"y":5,"g":8}

==========


プログラム自体は馬鹿馬鹿しいほど簡単なので,文字数を少なくするプログラムを目標にすると,以下のようなプログラムになった。

func = function(s) {
    x = matrix(unlist(strsplit(s, "")) == "G", byrow=TRUE, ncol=10)
    y = sapply(0:35, function(i) sum(x[i%%6+1:5, i%/%6+1:5]))
    i = which.max(y)-1
    cat(sprintf('{"x":%s,"y":%s,"g":%s}', i%/%6, i%%6, y[i+1]))
}
func("wwGwwwwwGGGwwwwwwwwwwwwwwwwwwwGwwwwwwwwwwwwwGwwGwwwGwwwwwwwwwwwGGwwwwwwwwwwwGwwwwwwwGGwwwwGwwwGGwGwG")

> func(s)
{"x":3,"y":5,"g":8}

コメント
この記事をはてなブックマークに追加

投影図から想像する立体

2016年11月01日 | ブログラミング

投影図から想像する立体
締め切りが 2016/11/01 10:00 AM なので,その 1 分後に投稿されるように予約

設問

3次元の物体を図面に表すとき、投影図を使うことがあります。
例えば、図の左側のようにブロックが積まれているとき、上面図・側面図・正面図を矢印の方向から見た図で考えると、
図の右側のようになります。

このとき、ブロックが見える位置を「1」、見えない位置を「0」とし、上の段から順に右のように表現します。

この上面図、側面図、正面図の表現が標準入力から与えられるとき、考えられるブロックの配置が何通りあるかを求め、
標準出力に出力してください。
なお、入力のサイズはいずれの面も3×3で固定とします。
また、ブロックの下には必ずブロックが必要で、空中に浮くブロックはないものとします。

例えば、上記の場合は他に以下のパターンがあり、全部で5通りですので、
以下のように出力します。

【入出力サンプル】
標準入力
[[1,1,1],[1,1,1],[1,0,1]]
[[0,0,1],[0,1,1],[1,1,1]]
[[1,0,0],[1,1,0],[1,1,1]]

標準出力
5

与えられた入力でブロックを配置できない場合は0を出力してください。

========

ブルートフォースでは 2 番目の入力例には無理だなあ

f = function(A, B, C) {
    B = B[3:1, 3:1]
    C = C[3:1, ]
    count = 0
    D = array(integer(27), dim=c(3,3,3))
    x = rbind(c(0L, 0L, 0L), c(1L, 0L, 0L), c(1L, 1L, 0L), c(1L, 1L, 1L))
    q = sum(A) == 27
    for (i1 in 1:4) {
        if (A[1,1] == 0 && i1 >= 2) break
        if (A[1,1] == 1 && i1 == 1) next
        D[1:3] = x[i1, ]
        if (sum(D[1:3] & !B[1:3])) break
        if (i1 == 1 && sum(B[1:3])) next
        for (i2 in 1:4) {
            if (A[2,1] == 0 && i2 >= 2) break
            if (A[2,1] == 1 && i2 == 1) next
            D[4:6] = x[i2, ]
            if (sum(D[4:6] & !B[4:6])) break
            for (i3 in 1:4) {
                if (A[3,1] == 0 && i3 >= 2) break
                if (A[3,1] == 1 && i3 == 1) next
                D[7:9] = x[i3, ]
                if (sum(D[7:9] & !B[7:9])) break
                for (i4 in 1:4) {
                    if (A[1,2] == 0 && i4 >= 2) break
                    if (A[1,2] == 1 && i4 == 1) next
                    D[10:12] = x[i4, ]
                    for (i5 in 1:4) {
                        if (A[2,2] == 0 && i5 >= 2) break
                        if (A[2,2] == 1 && i5 == 1) next
                        D[13:15] = x[i5, ]
                        for (i6 in 1:4) {
                            if (A[3,2] == 0 && i6 >= 2) break
                            if (A[3,2] == 1 && i6 == 1) next
                            D[16:18] = x[i6, ]
                            for (i7 in 1:4) {
                                if (A[1,3] == 0 && i7 >= 2) break
                                if (A[1,3] == 1 && i7 == 1) next
                                D[19:21] = x[i7, ]
                                for (i8 in 1:4) {
                                    if (A[2,3] == 0 && i8 >= 2) break
                                    if (A[2,3] == 1 && i8 == 1) next
                                    D[22:24] = x[i8, ]
                                    for (i9 in 1:4) {
                                        if (A[3,3] == 0 && i9 >= 2) break
                                        if (A[3,3] == 1 && i9 == 1) next
                                        D[25:27] = x[i9, ]
                                        if (all((apply(D, c(2, 3), sum) > 0) == A,
                                                (apply(D, c(1, 2), sum) > 0) == B,
                                                (apply(D, c(1, 3), sum) > 0) == C)) {
                                            count = count+1
                                        }
                                    }
                                }
                            }
                        }
                    }
                }
            }
        }
    }
    cat(count)
}

system.time({
f(matrix(c(1,1,1,1,1,1,1,0,1), byrow = TRUE, ncol = 3),  # 5
  matrix(c(0,0,1,0,1,1,1,1,1), byrow = TRUE, ncol = 3),
  matrix(c(1,0,0,1,1,0,1,1,1), byrow = TRUE, ncol = 3))
})
system.time({
f(matrix(c(1,1,1,1,1,1,1,1,1), byrow = TRUE, ncol = 3), # 4051
  matrix(c(1,1,1,1,1,1,1,1,1), byrow = TRUE, ncol = 3),
  matrix(c(1,1,1,1,1,1,1,1,1), byrow = TRUE, ncol = 3))
})
system.time({
f(matrix(c(1,0,1,0,1,0,1,0,1), byrow = TRUE, ncol = 3), # 17
  matrix(c(1,1,1,1,1,1,1,1,1), byrow = TRUE, ncol = 3),
  matrix(c(1,1,1,1,1,1,1,1,1), byrow = TRUE, ncol = 3))
})
system.time({
f(matrix(c(1,1,1,1,0,1,1,1,1), byrow = TRUE, ncol = 3), # 104
  matrix(c(0,0,0,1,1,1,1,1,1), byrow = TRUE, ncol = 3),
  matrix(c(0,0,0,1,1,1,1,1,1), byrow = TRUE, ncol = 3))
})
system.time({
f(matrix(c(1,1,1,0,1,0,1,1,1), byrow = TRUE, ncol = 3), # 0
  matrix(c(0,1,0,1,1,1,1,0,1), byrow = TRUE, ncol = 3),
  matrix(c(1,0,1,1,1,1,0,1,1), byrow = TRUE, ncol = 3))
})

コメント
この記事をはてなブックマークに追加

PPAP で遊々プログラミング

2016年10月30日 | ブログラミング

一つ目

I.HAVE.A = I.HAVE.AN = function(x) x
cat(rev(c(c(I.HAVE.A("PEN"), I.HAVE.AN("APPLE")),
rev(c(I.HAVE.A("PEN"), I.HAVE.A("PINEAPPLE"))))))

二つ目...わざわざ複雑そうに

PPAP = function(s, t, u, v) {
    cat1 = function(s, ...) cat(toupper(substr(s, 1, 1)), ...)
    cat2 = function(s) cat(sprintf("I have a%s %s.", ifelse(grepl(tolower(substr(s, 1, 1)), "aiueo"), "n", ""), s), "\n")
    cat3 = function(s) {cat2(s[1]); cat2(s[2]); cat(s[2:1], ".\n")}
    sapply(c(s, t, u, v, "\n"), cat1)
    sapply(list(c(s, t), c(u,v)), cat3)
    cat(sprintf("%2$s %1$s, %4$s %3$s.\n", s, t, u, v))
    cat(sprintf("%3$s %4$s, %2$s %1$s!\n", s, t, u, v))
}
PPAP("pen", "apple", "pen", "pineapple")
PPAP("ball", "baseball", "pen", "pencl")

コメント
この記事をはてなブックマークに追加

限られた数字で作れる数の、まんなかの値

2016年10月21日 | ブログラミング

締め切りが 10/21 10:00 AM というので,その 1 分後に投稿されるように予約

限られた数字で作れる数の、まんなかの値」を求める

【概要】
たとえば。
0, 4, 5以外の数を使わずに作ることができる3桁以下の非負の整数は、

    0, 4, 5, 40, 44, 45, 50, 54, 55, 400, 404, 405, 440, 444, 445, 450, 454, 455, 500, 504, 505, 540, 544, 545, 550, 554, 555

の27個です。中央は14番目の数なので、求める値は444になります。
このような感じで、桁数の最大値と使える数字を与えますので、作ることができる数の中央の値を求めてください。

【入出力】
入力は
3 045
こんな感じです。
空白の前が桁数の上限。空白の後が使える数字です。使える数字は重複なく昇順に並んでいます。

出力は、
444
のように、中央の値を出力してください。
作れる数が偶数個の場合には、中央に近い数を2つ、小さい順にコンマ区切りで並べてください。

【例】
入力     出力
3 045     444
1 012789     2,7
10 059     5555555555

【補足】
桁数の上限の最小値は 1 。最大値は 15 です
不正な入力に対処する必要はありません

================================

f = function(s) {
    g = function(m) {
        d = 0
        if (m > 0) {
            m = m-1 # 与えられる中央値は 1 から始まるので,j 進数で 0 から始めるために 1 を引く
            for (i in 1:99) {
                if (m == 0) break
                d[i] = x[m %% j + 1] # j 進数で各桁を決める。ただし,表記に使うのは x にある,実際に使える数字
                m = m %/% j
            }
        }
        paste(rev(d), collapse="")
    }
    options(scipen=100)
    s = unlist(strsplit(s, " "))
    k = as.numeric(s[1])
    x = as.numeric(unlist(strsplit(s[2], "")))
    if (length(x) == 1) { # 1 桁の場合は特別に扱う(簡単なのだけど)
        f = unique(as.numeric(sapply(1:k, function(y) paste(rep(x, y), collapse=""))))
        n = length(f)
        m = n %/% 2
        if (n %% 2 == 1) {
            cat(f[m+1])
        } else {
            cat(f[m+0:1], sep=",")
        }
    } else {
        j = length(x) # 使える数字の個数を j
        n = j^k # 作られる数の個数は 0 〜 (j^k)-1 の n 個。それぞれに 1 を足して 1 〜 n=j^k とする
        if (n %% 2 == 0) { # 真ん中の数は,要するに中央値。求め方は初歩の統計学の教科書に書いてある
            m = n/2
            cat(g(m), ",", g(m+1), sep="") # g() は十進数表記の中央値を j 進数で表記し直すだけ
        } else {
            m = (n+1) / 2
            cat(g(m))
        }
    }
}

> f("1 2") # 2
2
> f("2 2") # 2,22
2,22
> f("15 3") # 33333333
33333333
> f("1 0") # 0
0
> f("2 0") # 0
0
> f("15 04") # 44444444444444,400000000000000
44444444444444,400000000000000
> f("1 0246") # 2,4
2,4
> f("2 0246") # 26, 40
26,40
> f("15 0246") # 266666666666666,400000000000000
266666666666666,400000000000000
> f("15 12345789")
499999999999999,511111111111111
> f("15 123456789")
555555555555555
> f("15 0123456789") # 499999999999999,500000000000000
499999999999999,500000000000000

以下の 2 ケースの解が異なるのだが
■テストケース10 ・・・ 【入力値:15 12345789】【期待値:444444444444444,444444444444445】
■テストケース11 ・・・ 【入力値:15 123456789】【期待値:494949494949495】

コメント
この記事をはてなブックマークに追加

数列の判定

2016年10月14日 | ブログラミング

締め切りが 2016/10/24 10:00 AM なので,その 1 分後に投稿されるように予約。
当初の締め切りが延期されても,当方は一切関知しない。

【問題】
スペース区切りで 5個の10 進数が並んでいます。
与えられた数列が以下のいずれの数列(の、途中の連続した5項)なのかを判別するプログラムを書いてください。
記号     名称
G     等比数列
A     等差数列
F     フィボナッチ数
x     G, A, F のいずれにも該当しない

出力は、G, A, F, x のいずれかを出力してください。


【例】
入力     出力
1 2 4 8 16     G
1 2 3 4 5     A
3 5 8 13 21     F
1 2 123 1234 9999     x



【補足】
入力に含まれる値は、1以上、100億以下の整数です。
入力は、全て狭義単調増加列になっています。
不正な入力に対処する必要はありません。
G, A, F は大文字ですが、 x は小文字です。
1, 4, 5, 9, 14, 23, ... という数列はフィボナッチではありません。

f = function(s) {
    mul = function(x, y) { # 等比数列の各項が大きな値のとき,等比になっているかどうかは,ある項の二乗が前後の二項の積に等しいかどうかを見ないといけない
        x = rev(as.integer(unlist(strsplit(as.character(x), ""))))
        y = rev(as.integer(unlist(strsplit(as.character(y), ""))))
        l.x = length(x)
        l.y = length(y)
        z = integer(l.x+l.y)
        for (i in 1:l.x) {
            for (j in 1:l.y) {
                z[i+j-1] = z[i+j-1]+x[i]*y[j]
            }
        }
        for (i in 1:(l.x+l.y-1)) {
            z[i+1] = z[i+1] + z[i] %/% 10
            z[i] = z[i] %% 10
        }
        paste(rev(z), collapse="")
    }
    sqr = function(x) { # 二乗する関数
        mul(x, x)
    }
    d = as.numeric(unlist(strsplit(s, " ")))
    dif = diff(d)
    if (all(dif[1] == dif)) {
        cat("A")
    } else if (mul(d[1], d[3]) == sqr(d[2]) && # a[i-1] * a[i+1] == a[i]^2
               mul(d[2], d[4]) == sqr(d[3]) &&
               mul(d[3], d[5]) == sqr(d[4])) {
        cat("G")
    } else {
        n = 1:49
        fib = round((((1+sqrt(5))/2)^n - ((1-sqrt(5))/2)^n) / sqrt(5))
        i = which(fib == d[1])
        if (length(i) != 0 && all(d == fib[i[1]+0:4])) {
            cat("F")
        } else {
            cat("x")
        }
    }
}
f("160816114 194672138 235655746 285267482 345323794") # G
f("9845600625 9876856500 9908211600 9939666240 9971220736") # G
f("9715064832 9751864320 9788803200 9825882000 9863101250") # G
f("1234 6912 12590 18268 23946") # A
f("9999999990 9999999991 9999999992 9999999993 9999999994") # A
f("1 2000000000 3999999999 5999999998 7999999997") # A
f("13 21 34 55 89") # F
f("121393 196418 317811 514229 832040") # F
f("1134903170 1836311903 2971215073 4807526976 7778742049") # F
f("9983848527 9987006973 9990166419 9993326864 9996488308") # x
f("969323029 1568397607 2537720636 4106118243 6643838879") # x
f("1032569015 1670731762 2703300777 4374032539 7077333316") # x

コメント
この記事をはてなブックマークに追加

「トライアングル・メイズ」問題

2016年09月29日 | ブログラミング

締め切りが 2016/09/29 10:00 AM なので,その 1 分後に投稿されるように予約

設問

次のルールで生成される迷路を考えます。

レベル 1 の迷路から始めましょう。レベル 1 の迷路を M1 と呼びます。M1 は、3 つの点を L の字の形につないだものです。
左上の点がスタートで、右下の点がゴールです。

レベル 2 の迷路 M2 は、M1 を 3 つつなげて作ります。
M1 を 1 つ置き、さらにその上と右に 1 つずつ M1 をつなげたものです。
左上の点がスタートで、右下の点がゴールです。

レベル 3 の迷路 M3 は、M2 を 3 つつなげて作ります。
M2 を 1 つ置き、さらにその上と右に 1 つずつ M2 をつなげたものです。
左上の点がスタートで、右下の点がゴールです。

以降同様に、レベル k の迷路 Mk を、レベル k-1 の迷路 Mk-1 を 3 つつなげることで定義します。



さて、この迷路を、点と線をたどってスタートからゴールまで着く方法を考えます。
自然数 n に対し、迷路 Mn を最短距離でゴールするたどり方の数を F(n) と定義します。

例えば F(1) = 1,F(2) = 2 です。



同様に、F(3) = 6,F(4) = 42 です。
さらに、F(10) を 1000003(=10^6+3) で割った余りは 998593 となることが確かめられます。

標準入力から、自然数 n(1 ≦ n ≦ 10^8)が与えられます。
標準出力に F(n) を 1000003 で割った余り を出力するプログラムを書いてください。

=====

いつもの定石,すなわち,小さな問題の解を求め規則性を見つける。a(n) = a(n-1) * (a(n-1)+1)
ただし,真っ正直に計算するのでは n がばかでかいときには時間が掛かりすぎる。
ここでは,「1000003 で割った余り」を求めよということなので,答えの数列は一度ループにはまると,ループをぐるぐる回ることになるということ(これも定石)。
1~213 はループの外,214 項目の数値と同じ数値は1169 項目と同じになる(これは調べる)。
すなわち,214~1168, 1169~2123, ... は長さ955 のループ(同じ数字列の繰り返し)。

f = function(n) {
    m = 1168
    k = 1000003
    a = numeric(m)

    options(scipen=100)
    a[1] = f1 = 1
    a[2] = f2 = 2
    for (i in 3:m) {
        f3 = ((f2 %% k) * ((f2+1) %% k)) %% k
        a[i] = f3
        f1 = f2
        f2 = f3
    }
    if (n >= 214) {
        n = (n-214) %% 955 + 214
    }
    cat(a[n], " ")
}

> f(10)
998593  

> f(100000000)
342482  

コメント
この記事をはてなブックマークに追加

素数で作る天秤ばかり

2016年09月20日 | ブログラミング

素数で作る天秤ばかり

2016/09/20 10:00 AM なので,その 1 分後に投稿されるように予約

【問題】
天秤ばかりを使って重さを量りたいと考えています。
ただし、使えるおもりは重さが素数のものしかありませんでした。

m, n をともに正の整数とし、m 以下の素数すべてがおもりの重さとして1つずつ用意されているとき、n グラムの計り方が何通りあるかを求めてください。

例えば、m = 10, n = 2のとき、2, 3, 5, 7 のおもりが一つずつありますので、左右に以下のおもりを使った4通りがあります。
(量るものを左側に置いたとします。)

左側         右側
なし         「2」
「3」         「5」
「5」         「7」
「2」、「3」     「7」

標準入力から m と n がスペースで区切って与えられるとき、n グラムの計り方が何通りあるかを標準出力に出力してください。
(ただし、 m < 50 とします。)

【入出力サンプル】
標準入力
10 2

標準出力
4


==========

bin を毎回呼び出すとそのオーバーヘッドがあるので,結果をリスト bin.list に保存しておき,関数呼び出しをリストの参照に換えることで,実行時間は速くなる

bin = function(p) {
    retval = matrix(0L, nrow = 2^p, ncol = p)
    for (n in 1:p) {
        retval[, n] = rep(c(rep(0L, (2^p/2^n)), rep(1L, (2^p/2^n))), length = 2^p)
    }

    retval
}

f = function(m, n) {
    bin.list = sapply(1:14, bin)
    count = 0
    prime = c(2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47)
    prime = as.integer(prime[prime
    left.use = bin.list[[length(prime)]] # bin(length(prime))
    left = left.use %*% prime + n
    for (i in 1:(nrow(left.use) - 1)) {
        rest = prime[!left.use[i, ]]
        if (left[i]
            right.use = bin.list[[length(rest)]] # bin(length(rest))
            right = right.use %*% rest
            count = count + sum((left[i] == right))
        }
    }
    cat(count)
}


f(10, 2)  # 4
f(20, 10) # 90
f(30, 50) # 286
f(40, 30) # 3217
f(46, 3)  # 24946

コメント
この記事をはてなブックマークに追加

何回いえばわかるんだ(^_^)

2016年09月14日 | ブログラミング

> 第5回みえ県民意識調査は、各市町の選挙人名簿を使用した等間隔無作為抽出法により、標本を抽出しており、標本数10,000人に対して、有効回答数は5,236人でした。そのため、各属性において、実際の県全体と回答者の構成が異なる部分もあることから、以下にその概略をまとめています。

「等間隔無作為抽出法」というのは統計学用語にあるのかな?ふつうは「等間隔抽出法」とか「系統抽出法」じゃないか?

標本数10,000人」というのは大間違い。「調査対象が10,000人」ということ。「標本数」と「標本の大きさ(サンプルサイズ)」は大違い。

また,本文中で「回答数(サンプル数)」と何回も書いているが,「サンプル数」というのは統計用語にはない「標本の大きさ」とか「サンプルサイズ」というのがよい。

> 例えば、同じ調査を異なる調査対象で100回行った場合、95回以上の割合で同様の差が生じる場合は「統計的に有意な差がある」と表現し、90回以上の割合で同様の差が生じる場合は「統計的にある程度有意な差がある」と表現しています。

「統計学的に有意」の説明はいい加減。

> U>1.64の時、平均値の差は統計的に有意であると言える(危険率5%)

うおっ。何の断りもなく「片側検定」を行っているよ。

いい加減な用語を使っていると,内容の信頼性が低下する。

顧問として,「○○大学地域学部 教授 ○○○○」の名前があるんだが...

地元なんだから,三重大学の奥村晴彦先生に監修してもらうと良いのでは?

コメント
この記事をはてなブックマークに追加