オイラーのトーシェント関数

オイラーのトーシェント関数は、以下の形で表すことができます(正確な定義は オイラーのφ関数 - Wikipedia にあります)。

\varphi(n) = n (1 - \frac{1}{p_1}) (1 - \frac{1}{p_2}) \cdots (1 - \frac{1}{p_k})

この式の括弧内の分母を揃えます。

\varphi(n) = n (\frac{p_1 - 1}{p_1}) (\frac{p_2 - 1}{p_2}) \cdots (\frac{p_k - 1}{p_k})

プログラムで値を得るときは、この式を元に以下のように計算していました。

import std.experimental.all;

int totient(int n) {
    int ret = n;
    for (int p = 2; p * p <= n; p++) {
        if (n % p == 0) {
            ret = ret / p * (p - 1);
            while (n % p == 0) n /= p;
        }
    }
    if (n != 1) ret = ret / n * (n - 1);
    return ret;
}

void main() {
    writeln(iota(1, 20).map!totient);
    //=> [1, 1, 2, 2, 4, 2, 6, 4, 6, 4, 10, 4, 12, 6, 8, 8, 16, 6, 18]
}

npca2013年部誌_競技プログラミングと初等整数論入門(PDF) と読むと、 この関数を次のように計算していました。

int totient2(int n) {
    int ret = n;
    for (int p = 2; p * p <= n; p++) {
        if (n % p == 0) {
            ret -= ret / p; // ここが違う
            while (n % p == 0) n /= p;
        }
    }
    if (n != 1) ret -= ret / n; // ここが違う
    return ret;
}

こちらのほうが計算が単純になっていますね。この計算方法がどこから来るのか考えてみます。

最初の式に戻りますと、

\varphi(n) = n (1 - \frac{1}{p_1}) (1 - \frac{1}{p_2}) \cdots (1 - \frac{1}{p_k})

この式の最初の n(1 - \frac{1}{p_1}) を掛けると

\varphi(n) = (n - \frac{n}{p_1}) (1 - \frac{1}{p_2}) \cdots (1 - \frac{1}{p_k})

となります。最初の括弧内の項では、n から \frac{n}{p_1} を引いています。 次に最初の括弧内の項を A と置くと

\varphi(n) = A (1 - \frac{1}{p_2}) \cdots (1 - \frac{1}{p_k})

となり、A を展開すると

\varphi(n) = (A - \frac{A}{p_2}) \cdots (1 - \frac{1}{p_k})

となります。ここで行っているのは、A から \frac{A}{p_2} を引いています。
つまり、これは ret -= ret / p がやっていることなのですね。

参考

D言語: recurrence の練習

import std.experimental.all;

void main() {
    // 初項 1, 公差 2
    auto a = recurrence!((a, n) => a[n-1] + 2)(1);
    writeln(a.take(10));
    //=> [1, 3, 5, 7, 9, 11, 13, 15, 17, 19]

    // 初項 3, 公比 2
    auto b = recurrence!((a, n) => a[n-1] * 2)(3);
    writeln(b.take(10));
    //=> [3, 6, 12, 24, 48, 96, 192, 384, 768, 1536]

    // フィボナッチ数列
    auto fibs = recurrence!((a, n) => a[n-2] + a[n-1])(1, 1);
    writeln(fibs.take(10));
    //=> [1, 1, 2, 3, 5, 8, 13, 21, 34, 55]

    // 三角数
    auto t = recurrence!((a, n) => a[n-1] + n)(0).drop(1);
    writeln(t.take(10));
    //=> [1, 3, 6, 10, 15, 21, 28, 36, 45, 55]
}

参考

D言語: foreach の変数展開の仕組みが分からない

import std.experimental.all;

void main() {
    foreach (a; [10, 20, 30])
        writeln(a);
}

と書くと、配列の中身が foreach の変数 a に代入されて以下のように出力されます。

10
20
30

次に foreach の変数をひとつ追加します。 すると、その追加した変数には配列のインデックスが格納されます。

import std.experimental.all;

void main() {
    foreach (i, a; [10, 20, 30]) // 変数 i を追加した
        writeln(i, " ", a); // i はインデックス。a には配列の各要素
}

実行結果です。

0 10
1 20
2 30

さて、次は標準ライブラリの cartesianProduct 関数を使ったサンプルです。この関数は、引数で渡した配列の各要素の組み合わせを求める関数です。

import std.experimental.all;

void main() {
    foreach (a; cartesianProduct([1, 2], [3, 4]))
        writeln(a);
}

実行結果です。要素はタプルとしてまとめられます。

Tuple!(int, int)(1, 3)
Tuple!(int, int)(1, 4)
Tuple!(int, int)(2, 3)
Tuple!(int, int)(2, 4)

ここで foreach の変数をひとつ追加します。このとき、各変数にはタプルの各要素が代入されます。

import std.experimental.all;

void main() {
    foreach (a, b; cartesianProduct([1, 2], [3, 4]))
        writeln(a, " ", b);
}

実行結果です。

1 3
1 4
2 3
2 4

なるほど、各タプルの要素が展開されるのですね。 では、たとえば [tuple(1, 2), tuple(3, 4)] のような、タプルの配列があるとき foreach で各タプルの要素を変数で受けようとしても期待する結果にはなりません。

import std.experimental.all;

void main() {
    foreach (a, b; [tuple(1, 2), tuple(3, 4)])
        writeln(a, " ", b);
}

実行結果です。

0 Tuple!(int, int)(1, 2)
1 Tuple!(int, int)(3, 4)

上記のように、foeach の変数 a には配列のインデックスが格納されます。 配列と、cartesianProduct 関数のときで違いがでるのがよく分からないですね。

foreach の変数 a に 1 と 3 が代入されるようには書けないのですかね?

追記

twitter で教えていただきました。感謝です!

配列を zip でくるんで Input Range にするサンプル:

import std.experimental.all;

void main() {
    foreach (a, b; zip([tuple(1, 2), tuple(3, 4)])) {
        writeln(a, " ", b);
    }
}

実行結果です。

1 2
3 4

D言語: 2, 3, 5, 7... のシーケンスを作る

D 言語で、2, 3, 5, 7... のシーケンスを作る方法です。

素数をふるいにかけるとき

  • 2 で割り切れるか
  • 3, 5, 7... で割り切れるか

を試しますが、これらを分離せずにひとつのシーケンスとして扱えると便利ですね。

import std.experimental.all;

void main() {
    // 有限のサイズ
    auto a = chain([2], iota(3, 10, 2));
    writeln(a); // [2, 3, 5, 7, 9]

    // 無限のサイズ
    auto b = chain([2], recurrence!((a, n) => a[n-1] + 2)(3));
    writeln(b.take(10)); // [2, 3, 5, 7, 9, 11, 13, 15, 17, 19]
}

隣り合う要素でグループ分け

0 と 1 からなる配列があります。 たとえば [0, 1, 1, 0, 0, 1, 0, 0] の配列を [[0], [1, 1], [0, 0], [1], [0, 0]] のように 0 と 1 とでグループ分けする関数が必要になったので作ってみました。

def group(xs)
  return [] if xs.empty?
  g = []
  buf = [xs[0]]
  xs.drop(1).each {|e|
    same = yield buf[-1], e
    if same
      buf << e
    else
      g << buf
      buf = [e]
    end
  }
  g << buf
  g
end

実行結果です。

p group([1, 1, 0, 0, 0, 1, 0, 1, 1]) {|a, b| a == b }
#=> [[1, 1], [0, 0, 0], [1], [0], [1, 1]]

# 偶奇でグループ分け
p group([2, 4, 6, 1, 3, 6, 8]) {|a, b| a % 2 == b % 2 }
#=> [[2, 4, 6], [1, 3], [6, 8]]

Ruby で三角数の無限シーケンスを作る

まずは Enumerator で三角数を具体的に生成してみよう。

def triangles
  Enumerator.new do |y|
    sum = 0
    1.step {|i|
      sum += i
      y << sum
    }
  end
end

p triangles.take(10)
#=> [1, 3, 6, 10, 15, 21, 28, 36, 45, 55]

三角数を生成するところをもっと抽象化したい。
iterate を作ってみる。

def iterate(init)
  x = init
  Enumerator.new do |y|
    y << x
    loop { y << (x = yield x) }
  end
end

def triangles2
  x = 1
  iterate(1) {|sum| x += 1; sum + x }
end

p triangles2.take(10)
#=> [1, 3, 6, 10, 15, 21, 28, 36, 45, 55]

iterate の外側で変数(x)を定義しないようにできないかな。

def triangles3
  iterate([1, 2]) {|a, b| [a + b, b + 1] }.lazy.map {|e| e[0] }
end

p triangles3.take(10).to_a
#=> [1, 3, 6, 10, 15, 21, 28, 36, 45, 55]

iterate 内で三角数とインデックスを引き回すようにしてみた。
lazy の部分はきちんと理解していない。

Haskell には scanl という関数があるみたいだ。
Ruby で scanl を定義してみよう。

module Enumerable
  def scanl(z)
    Enumerator.new do |y|
      y << z
      x = z
      each {|e|
        y << (x += e)
      }
    end
  end
end

def triangles4
  2.step.scanl(1) {|a, b| a + b }
end

p triangles4.take(10)
#=> [1, 3, 6, 10, 15, 21, 28, 36, 45, 55]