総和 (Σ) 関数を書く

はじめに

数学は詳しくないのですが、「総和 というのがあって Σ 記号で表して...」というくらいは知っていて、「プログラミング言語で言うところの for 文に相当する」というのも理解していました。

• 総和 (Wikipedia)

そして調べ物をしていたら、Simula で書かれた関数 Sigma() というのが Wikipedia に載っていました。

• Simula (Wikipedia)

z= \sum_{i=1}^{100} \frac{1}{(i+a)^2}　

上記の式は

Z:= Sigma(i, 1, 100, 1 / (i + a) ** 2) ;

と書けるのだそうです。そしてその実装はどうなっているのかというと...

Real Procedure Sigma (l, m, n, u) ;
   Name l, u ;
   Integer l, m, n ;
   Real u ;
   Begin
   Real s ;
   l:= m ;
   While l <= n Do
      Begin
      s := s + u ;
      l := l + 1 ;
      End ;
   Sigma := s ;
   End ;

あー、簡単そうに見えて簡単に書けないやつだ Σ(ﾟдﾟlll)

See also:

• 名前呼び (Wikipedia)

実装

総和の計算をルーチンで書くなら難しくはないです。アルゴリズム解説の書籍では、大抵 for 文を使ったルーチンで書かれていると思います。実際、アルゴリズムの解説ならそれで用は足りますので。

program pgSigma;
{$APPTYPE CONSOLE}

uses
  System.Math;

function Sigma: Extended;
var
  i, lFrom, lTo: Integer;
  a: Extended;
begin
  lFrom := 1;
  lTo := 100;
  a := 2;
  result := 0.0;
  for i := lFrom to lTo do
    Result := Result + 1 / Power(i + a, 2);
end;

var
  z: Extended;
begin
  z := Sigma;
  Writeln(z: 1 : 2);
end.

Delphi (Pascal) にはべき乗の演算子がないので Power() 関数を使っています。

でも、Simula の Sigma() 関数と同等な、

• 式を関数のパラメータとして与えることができる
• 副作用を伴わない

上記を満たす関数を書こうとすると、途端に頭を抱える事になります。「そんなの簡単じゃないか？」と思った方は次のような実装を考えたのだと思います。

function Sigma(var i: Integer; const aFrom, aTo: Integer; aExp: Extended): Extended;
begin
  Result := 0.0;
  i := aFrom;
  while i <= aTo do // 言語によっては for にローカル変数しか使えないので while で
  begin
    Result := Result + aExp;
    Inc(i); // i++
  end;
end;

次のような式の場合、

z=\sum_{i=1}^{5} 2 = 2 + 2 + 2 + 2 + 2

呼び出しはこんな感じで。

  var i: Integer; 
  var z := Sigma(i, 1, 5, 2);
  Writeln(z: 1 : 2);

確かに結果は正しく 10 となりました。

10.00

でも、次のような式の場合、

z=\sum_{i=1}^{5} i = 1 + 2 + 3 + 4 + 5

  var i: Integer; 
  var z := Sigma(i, 1, 5, i);
  Writeln(z: 1 : 2);

結果は正しくなりません。

0.00

何故なら上記コードは次のコードと同等だからです (変数 i が 0 で初期化される場合)。

  var i: Integer; 
  var z := Sigma(i, 1, 5, 0);
  Writeln(z: 1 : 2);

では 4 番目のパラメータを変数パラメータ (参照渡し) にすると、

function Sigma(var i: Integer; const aFrom, aTo: Integer; var e: Integer): Extended;
begin
  Result := 0.0;
  i := aFrom;
  while i <= aTo do
  begin
    Result := Result + e;
    Inc(i);
  end;
end;

結果は正しくなりましたが...

15.00

当然、式が使えなくなります。

z=\sum_{i=1}^{5} (i+1) = 2 + 3 + 4 + 5 + 6

  var i: Integer; 
  var z := Sigma(i, 1, 5, i + 1); // <- 式が使えない
  Writeln(z: 1 : 2);

これが実装が簡単でない理由です。

僕が用意した解答を見ずに、自分でも考えてみたいという方は、この辺でブラウザのスクロールを止めてください。

See also:

• 11.1.1. パラメータリスト (パラメータ並び / Parameter list) (Qiita)
• 11.2.2. 副作用 (Side effect) (Qiita)
• 12.3. 手続き Write と Writeln

・Delphi での実装

こんな感じのユニットを作ってみました。

uSigma.pas

unit uSigma;

interface

type
  TSigmaExp = reference to function (i: Integer): Extended;

function Sigma(var i: Integer; const aFrom, aTo: Integer; aExp: TSigmaExp): Extended;

implementation

function Sigma(var i: Integer; const aFrom, aTo: Integer; aExp: TSigmaExp): Extended;
begin
  Result := 0.0;
  i := aFrom;
  while i <= aTo do
  begin
    Result := Result + aExp(i);
    Inc(i);
  end;
end;
end.

メソッド参照型を使って、計算式はそこに書くようにしました。こうしておけば Sigma() に無名メソッドが渡せます。

例1: Z:= Sigma(i, 1, 100, 1 / (i + a) ** 2) ;

こうかな？

z= \sum_{i=1}^{100} \frac{1}{(i+a)^2}　

program pgSigma;
{$APPTYPE CONSOLE}
uses
  System.Math, uSigma;

begin
  var i: Integer;
  var a := 2; // 任意の数
  var z := Sigma(i, 1, 100,
    function (i: Integer): Extended
    begin
      result := 1 / Power(i + a, 2);
    end);
  Writeln(z: 1 : 2);
end.

a=2 の時の結果は次のようになりました。

0.39

例2: Z:= Sigma(i, 1, 5, i + 1) ;

z=\sum_{i=1}^{5} (i+1)　

pgSigma.dpr

program pgSigma;
{$APPTYPE CONSOLE}
uses
  uSigma;

begin
  var i: Integer;
  var z := Sigma(i, 1, 5,
    function (i: Integer): Extended
    begin
      result := i + 1;
    end);
  Writeln(z: 1 : 2);
end.

こちらは検算が簡単ですね (w

20.00

ちゃんと 15.00 が表示されました。

・標準 Pascal での実装

つまりは名前渡しができなくとも、関数のパラメータとして関数が渡せれば総和関数が書けます (コールバック関数など)。パラメータとして無名関数が渡せればもっとスマートになります。

標準 Pascal も関数パラメータに対応しているため、総和関数が書けます。

pgSigma.pas

program pgSigma(Output);

var
  i: Integer;
  z, a: Real;

  function Sigma(var i: Integer; aFrom, aTo: Integer; function f(i: Integer): Real): Real;
  var
    v: Real;
  begin
    v := 0.0;
    i := aFrom;
    while i <= aTo do
    begin
      v := v + f(i);
      i := i + 1;
    end;
    Sigma := v;
  end; { Sigma }

  function SigmaExp1(i: Integer): Real;
  begin
    SigmaExp1 := 1 / Exp(Ln(i + a) * 2);
  end; { SigmaExp1 }

  function SigmaExp2(i: Integer): Real;
  begin
    SigmaExp2 := i + 1;
  end; { SigmaExp2 }

begin
  { Sigma #1 }
  a := 2;
  z := Sigma(i, 1, 100, SigmaExp1);
  Writeln(z : 1 : 2);
  { Sigma #2 }
  z := Sigma(i, 1, 5, SigmaExp2);
  Writeln(z : 1 : 2);
end.

標準 Pascal には Power() 関数もないため、べき乗相当を Exp() と Ln() を使って計算しています。ただし、Xⁿ の時、X > 0 である必要があります。古い Delphi の Standard 版も Math ユニットが含まれないため、Power() 関数が使えませんでした。

式を分離して関数にしなくてはならないのなら、総和関数を 2 度書いても大して手間は変わらない気がするので、個人的にはイマイチに思えます。

See also:

• 11.2.1. 関数パラメータ (Functional parameters)
• System.Math (DocWiki)
• Pascal (J&W) - Wirth 先生の邦訳本を読んでみる (Qiita) ¹
• #12 コンピュータアルゴリズム事典 - 技術評論社の Software Technology シリーズから Pascal 関連の書籍を読んでみる (Qiita) ²

追記

後で『J&W (第4版)』の「11.2.1. 関数パラメータ」を確認したら、そこの例が 「プログラム 11.9 級数の連続和の表を書く」 で、そのコードの中に Sigma() 関数がありました！当然ですが、内容もほぼ同じでした。

function Sigma(function F(X: Real): Real; Lower, Upper: Integer): Real;
var
  Index: Integer;
  Sum: Real;
begin
  Sum := 0.0;
  for Index := Lower to Upper do
    Sum := Sum + F(Index);
  Sigma := Sum;
end; { Sigma } 

／(^o^)＼ﾅﾝﾃｺｯﾀｲ

「11.2.1. 関数パラメータ (原著では「11.B.1. Functional parameters」)」のサンプルコード (プログラム 11.9) は『J&W (初版・第2版)』と『J&W (第3版・第4版)』で内容が異なっているのですが、第2版の方しか頭に入ってなかったみたいです。

おわりに

今回記述した総和関数ですべての Σ 記号のある式を計算できる訳ではありませんが、Simula の Sigma() 関数と同等のものは書けたように思います (何か勘違いしてなきゃいいけど...)。もっとスマートに書ける方法をご存じな方がいらっしゃいましたら、ぜひ教えてください m(_ _)m

あと、他の言語では総和関数をどう書けるのか知らないので (標準で持っていたり？) 、他の言語の状況もコメ欄で教えて頂けると幸いです。

See also:

• Delphi 用に FizzBuzz コレクションを作る (Qiita)
• ホクホクのイモを Delphi で書いてみる (Qiita)
• Delphi で湯婆婆を実装してみる (Qiita)
• Delphi で前前前世 (Qiita)

1. 標準 Pascal で書いた Power() 関数が載っています (プログラム 11.8)。 ↩

2. 標準 Pascal で書いた Power() 関数が載っています。(累乗 <1>・<2>・<3>) ↩