DouKaku?beta

challenge 正整数のゲーデル数化？ 

正の整数 n を引数としてとり， 2^d1 * 3^d2 * 5^d3 ... * pk^dk を返す関数
goedel を定義してください．

ただし，n を10進表現で k 桁の数としたときの各桁の数が数列 [d1,d2,d3,...,dk]
をなすとし，dk が 1 の位，d1 が 10^(k-1) の位です．また，pk は k番目の素数です．

goedel   9  ⇒ 2^9             ⇒  512
goedel  81  ⇒ 2^8 * 3^1       ⇒  768
goedel 230  ⇒ 2^2 * 3^3 * 5^0 ⇒  108

see: ゲーデル数

post comment

Posted feedbacks - Nested

Flatten Hidden 
全力で問題を勘違いしていた。
素数を生成する部分を追加してみた。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
def CreatePrimes(n):
    primes = []
    c = 2
    while len(primes) < n:
        primes.append(c)
        for x in primes[:-1]:
            if c%x == 0:
                primes.pop();
                break
        c += 1
    return primes

def goedel(n):
    result = 1
    primes = CreatePrimes(len(str(n)))
    
    for x in str(n):
        result *= primes[0]**int(x)
        primes = primes[1:]
    return result

print goedel(9)
print goedel(81)
print goedel(230)
print goedel(12345)
こんな感じでしょうか。IEnumerator<int> がなんか気持ち悪いですが…… 
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
using System;
using System.Collections.Generic;
static class Program {
    static void Main() {
        Console.WriteLine(Goedel(9));   // 512
        Console.WriteLine(Goedel(81));  // 768
        Console.WriteLine(Goedel(230)); // 108
    }
    static double Goedel(int num) {
        double goedel = 1;
        IEnumerator<int> p = Prime();
        foreach(char c in num.ToString()) {
            p.MoveNext();
            goedel *= Math.Pow(p.Current, c - '0');
        }
        return goedel;
    }
    static IEnumerator<int> Prime() {
        yield return 2;
        List<int> primes = new List<int>();
        for(int num = 3; ; num += 2) {
            foreach(int p in primes) {
                if (num % p == 0)
                    goto SKIP;
            }
            primes.Add(num);
            yield return num;
        SKIP:;
        }
    }
}
n を 10 桁までに限定しています。 
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
#include <stdlib.h>
#include <stdio.h>
#include <math.h>

int goedel( int n ){
   static int prime[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29};
   char buf[11];
   char *dk = buf;
   int *pk = prime;
   int result = 1;
   snprintf(dk,11,"%d",n);
   while(*dk){
      result *= pow( *pk++, *dk++ -'0' );
   }
   return result;
}

void main ( void ){
   printf("%d\n", goedel(9));
   printf("%d\n", goedel(81));
   printf("%d\n", goedel(230));
}

もうちょっと簡単に出来そう

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
use strict;
use warnings;
use Math::Prime::XS qw(is_prime);
use List::Util qw(reduce);
use List::MoreUtils qw(pairwise);

sub primes {
    my $number = shift;
    my @primes = (2);
    my $i      = 3;
    while ( @primes < $number ) {
        push( @primes, $i ) if is_prime($i);
        $i++;
    }
    return @primes;
}

sub goedel {
    my $number = shift;
    my @num    = split //, $number;
    my @primes = primes( scalar(@num) );
    reduce { $a * $b } pairwise { $a**$b } @primes, @num;
}

print goedel(9);
print goedel(81);
print goedel(230);
Squeak Smalltalk で。

指定した数未満の素数列を生成する、組み込みの Integer class>>#primesUpTo: を流用して富豪的に。 
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
| primesOfSize goedel |

primesOfSize := [:n |
    | m primes  |
    m := 0.
    [(primes := Integer primesUpTo: (10 raisedTo: (m := m + 1))) size >= n] whileFalse.
    primes first: n].

goedel := [:num |
    | digits size |
    digits := (num asString as: Array) collect: [:char | char asString asInteger].
    size := digits size.
    primes := primesOfSize value: size.
    (1 to: size) inject: 1 into: [:goe :idx | goe * ((primes at: idx) raisedTo: (digits at: idx))]].

goedel value: 9.   "=> 512 "
goedel value: 81.   "=> 768 "
goedel value: 230.   "=> 108 "

やっぱり一度文字列化するの一番効率的で簡単みたいです。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
import Data.Char

primes = sieve [2..] 
  where 
  sieve (x:xs) = x : sieve [y| y<-xs, y `mod` x /= 0]

goedel x = product $ zipWith (^) primes (map digitToInt (show x))

{-
*Main> goedel 9
512
*Main> goedel 81
768
*Main> goedel 230
108
-}

10桁限定。それ以上は整数がoverflowします。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
fun prime n a 0 = a
  | prime n a i =
  if n mod 2 <> 0 andalso n mod 3 <> 0 then
    prime (n + 1) (a @ [n]) (i - 1)
  else
    prime (n + 1) a i


fun goedel n =
let
  val p = prime 5 [2, 3] 10
  val a = map (valOf o Int.fromString o str) ((explode o Int.toString) n)
in
  floor (ListPair.foldl (fn (x, y, z) => Math.pow (real x, real y) * z) 1.0 (p, a))
end

val printInt = println o Int.toString;

printInt (goedel 9);
printInt (goedel 81);
printInt (goedel 230)

素数生成が間違っていたので修正。 ついでに多倍長整数を使って、Overflowしないようにしてみた。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
fun prime n =
let
  val a = List.tabulate (n - 1, fn x => x + 2)

  fun loop [] = []
    | loop (x::xs) =
      x :: loop (List.filter (fn i => i mod x <> 0) xs)
in
  loop a
end


fun goedel (n : IntInf.int) =
let
  open IntInf

  val p = map fromInt (prime 100)
  val a = map (valOf o Int.fromString o str) ((explode o toString) n)
in
  toString (ListPair.foldl (fn (x, y, z) => pow (x, y) * z) 1 (p, a))
end


fun println s = print (s ^ "\n")

val _ = println (goedel 9)
val _ = println (goedel 81)
val _ = println (goedel 230)
val _ = println (goedel 19425463134)

あれ、これだと

goedel 230  ⇒ 2^2 * 3^3 * 5^0 ⇒  108
goedel 23   ⇒ 2^2 * 3^3 ⇒  108

なので、ゲーデル数と元の数との一対一対応が取れないのでは?

こうすればOKだけど

goedel' 230  ⇒ 2^(2+1) * 3^(3+1) * 5^(0+1) ⇒  3240
goedel' 23   ⇒ 2^(2+1) * 3^(3+1) ⇒  648

Dan the Goedel Numberer
ゲーテル数の定義とは関係ないですが、
一対一対応を素直に取れるように定義を変更するなら、
逆向きに並べるのが単純でいい気がします。

goedel 230  ⇒ 2^2 * 3^3 * 5^0 ⇒  108
goedel  23  ⇒ 2^2 * 3^3 ⇒  108

ではなく、

goedel 230  ⇒ 2^0 * 3^3 * 5^2 ⇒  675
goedel  23  ⇒ 2^3 * 3^2 * 5^0 ⇒  72

これで、位取りの意味的に上位の0は無視することと、
0乗した値が無視できることが対応づけられます。

goedel_number()の二番目の引数にnonzeroを指定することで#4657で指摘した「仕様バグ」にも対応しています。

Dan the Goedel Numberer 

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
#!/usr/local/bin/perl
use strict;
use warnings;

my @primes = ( 2, 3 );

sub fill_primes {
    my $nprimes = shift;
  ODD:
    for ( my $n = $primes[ @primes - 1 ] + 2 ; @primes < $nprimes ; $n += 2 ) {
        for my $d (@primes) {
            last if $d * $d > $n;
            next ODD unless $n % $d;
        }
        push @primes, $n;
    }
}

sub goedel_number {
    my $n      = shift;
    my $offset = shift || 0;
    my @d      = ( $n =~ /(\d)/g );
    fill_primes( scalar @d );
    use bigint;
    my $result = 1;
    $result *= $primes[$_]**( $d[$_] + $offset ) for ( 0 .. @d - 1 );
    $result;
}

local $\ = "\n";
local $, = ", ";
print goedel_number( $ARGV[0] || 1234567890, $ARGV[1] );

#4658の素直な移植。

Dan the Goedel Numberer 

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
var primes = [2,3];
function fill_primes(nprimes){
  function f(n){
    for (var p = 0, l = primes.length; p < 0; p++){
      if (p*p > n)    return;
      if (n % d == 0) break;
    }
    primes[primes.length] = n;
  }
  for (var n = primes[primes.length-1]+2; primes.length < nprimes; n += 2)
    f(n);
}
function goedel_number(n, offset){
  if (! offset) offset = 0;
  var d = [];
  n.toString().replace(/\d/g, function(m){
    d[d.length] = parseInt(m);
  });
  fill_primes(d.length);
  var result = 1;
  for (var i = 0, l = d.length; i < l ; i++) {
    result *= Math.pow(primes[i], d[i]+offset);
  }
  return result;
}
/*
p(goedel_number(23));
p(goedel_number(230));
p(goedel_number(23,1));
p(goedel_number(230,1));
*/

車輪は発明しない方向で。

1
2
3
4
5
require "mathn"
def goedel(n)
  prime=Prime.new
  n.to_s.split(//).inject(1){|r,k| r*=prime.next**k.to_i}
end

素朴に。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
import itertools

def is_prime(n):
    for i in itertools.count(2):
        if i * i > n:
            return True
        if n % i == 0:
            return False

def primes():
    for n in itertools.count(2):
        if is_prime(n):
            yield n

def numbers(n):
    a = []
    while n:
        a.append(n % 10)
        n /= 10
    return reversed(a)

def goedel(n):
    result = 1
    for pk, dk in itertools.izip(primes(), numbers(n)):
        result *= (pk ** dk)
    return result

def main():
    for n in (9, 81, 230):
        print goedel(n)

if __name__ == '__main__':
    main()

素数を作るところをもっとなんとかしたい・・・。 とりあえず出来たのでUP。

see: エラトステネスの篩 - Wikipedia 

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
<?php
function prime($len)
{
    $result = array(2);
    for($n = 2; $n <= $len; $n++) {
        if($n % 2 == 0) {
            continue;
        }
        for($i = 3 ; pow($i, 2) <= $n ; $i += 2) {
            if($n % $i == 0)
                continue 2;
        }
        array_push($result, $n);
    }
    return $result;
}

function goedel($n)
{
    $n = strval($n);
    $result = 1;
    $prime = prime($n);
    for($i = 0; $i < strlen($n); $i++) {
        $result *= pow($prime[$i], $n[$i]);
    }
    return $result;
}

特にひねりのない実装。
素数はエラトステネスのふるいを使ったジェネレータにしてますが、
工夫のない作りなのでたぶん指数オーダーです。

最近は関数型至上主義も薄れてきて、リファレンスや for 文なども
わりと平気で使えるようになってきました。
外に影響出さないなら、中で何やってたって良いんです :-p


% ocaml goedel.ml
goedel 9 => 2^9 => 512
goedel 81 => 2^8 * 3^1 => 768
goedel 230 => 2^2 * 3^3 * 5^0 => 108

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
let gen_prime () =
   let now = ref 2 in
   let buf = ref [2] in
   let rec is_prime = function
      | 2 -> true
      | n ->
           if List.exists (fun i -> n mod i = 0) !buf
           then false
           else (buf := n::!buf; true)
   in
   let rec next_prime i =
      let n = !now in
      incr now;
      if is_prime n
      then Some n
      else next_prime i
   in
   Stream.from next_prime

let goedel n =
   let pow n = function
      | 0 -> 1
      | m ->
           let rec loop acc = function
              | 1  -> acc
              | m' -> loop (acc * n) (m' - 1)
           in
           loop n m
   in
   let s      = string_of_int n in
   let buf    = Buffer.create 10 in
   let add    = Buffer.add_string buf in
   let result = ref 1 in
   let primes = gen_prime () in
   for i = 0 to (String.length s) - 1 do
      let m = int_of_string (String.make 1 s.[i]) in
      let n = Stream.next primes in
      if Buffer.length buf > 0 then add " * ";
      List.iter add [string_of_int n; "^"; string_of_int m];
      result := !result * (pow n m)
   done;
   (Buffer.contents buf, !result)

let test () =
   let samples = [9; 81; 230] in
   List.iter begin fun n ->
      let s, r = goedel n in
      Printf.printf "goedel %d => %s => %d\n" n s r
   end samples

let () = if not !Sys.interactive then test ()
なんとなく大きな数対応版 (Num モジュール使用)。
素数の方はたかだか数十番目くらいしか必要にならないので int のまま。

% ocaml nums.cma goedel.ml
goedel 9 => 2^9 => 512
goedel 81 => 2^8 * 3^1 => 768
goedel 230 => 2^2 * 3^3 * 5^0 => 108
goedel 1256 => 2^1 * 3^2 * 5^5 * 7^6 => 6617756250
goedel 34721 => 2^3 * 3^4 * 5^7 * 7^2 * 11^1 => 27286875000
goedel 542673 => 2^5 * 3^4 * 5^2 * 7^6 * 11^7 * 13^3 => 326393949142783922400
goedel 19425463134 => 2^1 * 3^9 * 5^4 * 7^2 * 11^5 * 13^4 * 17^6 * 19^3 * 23^1 * 29^3 * 31^4 => 475616767259341262526341725017954602561250 
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
open Num

let gen_prime () =
   let now = ref 3 in
   let buf = ref [2] in
   let rec is_prime = function
      | 2 -> true
      | n ->
           let limit = n / 2 in
           if List.exists begin fun i ->
              if i > limit
              then false (* skip *)
              else n mod i = 0
           end !buf
           then false
           else (buf := n::!buf; true)
   in
   let rec next_prime = function
      | 0 -> Some 2
      | i ->
         let n = !now in
         now := !now + 2;
         if is_prime n
         then Some n
         else next_prime i
   in
   Stream.from next_prime

let goedel num =
   let s      = string_of_num num in
   let buf    = Buffer.create 10 in
   let add    = Buffer.add_string buf in
   let result = ref (Int 1) in
   let primes = gen_prime () in
   for i = 0 to (String.length s) - 1 do
      let cs = String.make 1 s.[i] in
      let m  = num_of_string cs in
      let n  = Stream.next primes in
      if Buffer.length buf > 0 then add " * ";
      List.iter add [string_of_int n; "^"; cs];
      result := !result */ (power_num (Int n) m)
   done;
   (Buffer.contents buf, !result)

let test () =
   let samples = [9; 81; 230; 1256; 34721; 542673; 19425463134] in
   let pr ch n = output_string ch (string_of_num n) in
   List.iter begin fun n ->
      let s, r = goedel (Int n) in
      Printf.printf "goedel %a => %s => %a\n" pr (Int n) s pr r
   end samples

let () = if not !Sys.interactive then test ()

とりあえず，「仕様バグ」は無視して題意どおりの実装。

$ java Goedel 230
108

$ java Goedel 123456789
52713275010243038997421106186697438702252144407250
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
import java.math.BigInteger;
import java.util.Vector;

public class Goedel {
    private static final int CERTAINTY = 24;
    public static void main(String[] args) throws Exception {
        Vector<BigInteger> primes = new Vector<BigInteger>();
        for (int i = 2; primes.size() < args[0].length(); i++) {
            BigInteger bi = BigInteger.valueOf(i);
            if (bi.isProbablePrime(CERTAINTY)) primes.add(bi);
        }
        BigInteger result = BigInteger.ONE;
        for (int i = 0; i < args[0].length(); i++)
            result = result.multiply(primes.get(i).pow(Integer.valueOf(args[0].substring(i,i+1))));
        System.out.println(result);
    }
}

なんとなくstreamを使ってみた。 sieveとprimesはSICPのやつ。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
(use util.stream)
(use gauche.collection)

(define (integers-starting-from n)
  (stream-cons n (integers-starting-from (+ n 1))))

(define (sieve stream)
  (define (divisible? x y) (= (remainder x y) 0))
  (stream-cons
   (stream-car stream)
   (sieve (stream-filter
           (lambda (x) (not (divisible? x (stream-car stream))))
           (stream-cdr stream)))))

(define primes (sieve (integers-starting-from 2)))

(define (number->digit-stream n)
  (stream-map
   (lambda (c) (- (char->integer c) (char->integer #\0)))
   (string->stream (number->string n))))

(define (goedel n)
  (apply * (stream->list
            (stream-map (lambda (pk nk) (expt pk nk))
                        primes (number->digit-stream n)))))

所詮は浮動小数点数なので大きなnについての精度は無考慮。

primes(100)は100までの素数25個を返す。 25と決めうちしているのは気持ち悪いが、double型が正確に表わせるのは高々十数桁であり、その範囲で使うぶんにはとりあえずはOKとする。 桁数に応じた数の素数をとりだしたい場合は、第4行を次のように置き換えればよい:

if n < 6
  p = primes(11);  % Five smallest primes
else
  u = n*(log(n) + log(log(n)));  % An upper bound of the value of the nth prime
  p = primes(u);
end

上界uは ピエール・デザルトの定理 による。

1
2
3
4
5
function g = goedel(n)
  s = num2str(n) - '0';
  n = length(s);
  p = primes(100);
  g = prod(p(1:n).^s);

で、結局ゲーデル数って何ですか?(w

see: 貧脚レーサーのサボり日記 

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
//http://ja.doukaku.org/100/　投稿用
using System;
using System.Collections.Generic;
class Program {
　　static void Main(string[] args) {
　　　　Console.WriteLine(goedel(9));
　　　　Console.WriteLine(goedel(81));
　　　　Console.WriteLine(goedel(230));
　　　　Console.ReadLine();
　　}

　　static double goedel(int n) {
　　　　string tmpStr = n.ToString();
　　　　double r = 1;
　　　　List<double> prime = new List<double>(new double[] { 2 });
　　　　for(int i = 3; ; i += 2) {
　　　　　　bool isPrime = true;
　　　　　　for(int j = 3; j < i; j++) {
　　　　　　　　if(i % j == 0) {
　　　　　　　　　　isPrime = false;
　　　　　　　　　　break;
　　　　　　　　}
　　　　　　}
　　　　　　if(isPrime) prime.Add(i);
　　　　　　if(prime.Count >= tmpStr.Length) break;
　　　　}
　　　　for(int i = 0; i < tmpStr.Length; i++) {
　　　　　　r *= Math.Pow(prime[i], (double.Parse(tmpStr[i].ToString())));
　　　　}
　　　　return r;
　　}
}

とりあえずお題通りに（桁が変わると一意じゃないのを見越してタイトルに？が付いてるのと解釈）。 とりあえず解が出せるのは5桁までかな....6桁だと64bitでも溢れる気がする。 桁数が少ないので素数は計算せず。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
#include <complex>

static unsigned int prime[] = {2, 3, 5, 7, 11};

unsigned int digit(unsigned int n, unsigned int base)
{
    unsigned d = 1;
    for (; n > (base - 1); n /= base) { d++; }
    return d;
}

unsigned long long goedel(unsigned int n)
{
    unsigned long long ans = 1;
    for (unsinged int d = digit(n, 10); d > 0; d--, n /= 10)
        ans *= (unsigned long long)std::pow((double)prime[d - 1], (int)(n % 10));
    return ans;
}

計算違い。4桁で既に溢れてる....orz