多忙な日々は続く
2005年5月27日(金) 15:00多忙な日々が続いていますので、 日記更新が滞っています。 でも元気ですのでご心配なく。 お祈りください。
結城浩の日記
多忙な日々が続いていますので、 日記更新が滞っています。 でも元気ですのでご心配なく。 お祈りください。
神さま、 今日という日を与えてくださってありがとうございます。 いままでになく、これからもない、 最高の一日を与えてくださることを信じ、感謝します。 あなたが与えてくださっている人生を感謝します。 時間を、家族を、友人を感謝します。 神さま、あなたはなんとすばらしい方でしょう。 あなたは、必要なものをすべて備え、与えてくださる方です。 しかもそれを最善の時に与えてくださいます。感謝します。 あなたは命を与え、救いを与えてくださいました。 あなたの愛を感謝します。 主よ、あなたが与えてくださる恵みにもっともっと感謝することができるように、 物事をよく見る力を与えてください。 また、自分が語るだけではなく、あなたの言葉、 そしてほかの人の言葉に耳を傾けることができますように。
十字架に感謝し、復活に感謝し、 イエスさまのお名前で祈ります。
アーメン!
元気ですが、何かと多忙で、日記頻度が落ちております。
以下の文章は、LaTeXのソースになっています。
\documentclass[fleqn]{jsarticle}
\title{三角関数の離散系バージョン}
\author{結城浩\footnote{http://www.hyuki.com/\ Hiroshi Yuki \copyright\ 2005, All rights reserved.}}
\date{2005年5月}
\usepackage{amsmath}
\usepackage{euler}
\begin{document}
% \thispagestyle{empty}
\maketitle
% \tableofcontents
\par
この文章は、
「離散系バージョンの関数探し」\footnote{http://www.hyuki.com/story/diffsum2.html}の続きです。
「離散系バージョンの関数探し」では、
三角関数の離散系バージョンはあまり納得が行かない結果で終わりました。
そこで、級数展開を使って求められないかどうかを調べてみます。
\par
まず、
$\cos x$と$\sin x$の級数展開は次の通りです。
\begin{eqnarray*}
\cos x = + \frac{x^0}{0!} - \frac{x^2}{2!}
+ \frac{x^4}{4!} - \frac{x^8}{8!} + \cdots
= \sum_{k=0}^{\infty} (-1)^k \frac{x^{2k+0}}{(2k + 0)!}\\
\sin x = + \frac{x^1}{1!} - \frac{x^3}{3!}
+ \frac{x^5}{5!} - \frac{x^7}{7!} + \cdots
= \sum_{k=0}^{\infty} (-1)^k \frac{x^{2k+1}}{(2k + 1)!}\\
\end{eqnarray*}
\par
$x^n$の部分を機械的に下方階乗冪$x^{\underline{n}}$に書き換えると、次のようになります。
\begin{eqnarray*}
c(x) = + \frac{x^{\underline{0}}}{0!} - \frac{x^{\underline{2}}}{2!}
+ \frac{x^{\underline{4}}}{4!} - \frac{x^{\underline{8}}}{8!} + \cdots
= \sum_{k=0}^{\infty} (-1)^k \frac{x^{\underline{2k+0}}}{(2k + 0)!}\\
s(x) = + \frac{x^{\underline{1}}}{1!} - \frac{x^{\underline{3}}}{3!}
+ \frac{x^{\underline{5}}}{5!} - \frac{x^{\underline{7}}}{7!} + \cdots
= \sum_{k=0}^{\infty} (-1)^k \frac{x^{\underline{2k+1}}}{(2k + 1)!}\\
\end{eqnarray*}
\par
項別に差分を取って確かめると、
$c(x)$と$s(x)$は以下の差分方程式を満たします(収束するかどうかの吟味は無視しています)。
\begin{eqnarray*}
\Delta c(x) & = & - s(x) \\
\Delta s(x) & = & c(x) \\
\end{eqnarray*}
うん、これなら三角関数の離散系バージョンがきちんと出るかもしれませんね。
どんな関数なんでしょう。わくわく。
振る舞いを把握するために、具体的な値を計算してみましょう。
はじめの20項までを計算して、$x = 0, 1, 2, ..., 15$のときの$c(x)$と$s(x)$を求めるPerlのプログラムを書きました。
\begin{verbatim}
use strict;
use bignum;
my $maxn = 20;
for my $n (0..15) {
my $c = &c($n);
my $s = &s($n);
printf("%3d %8g %8g\n", $n, $c, $s);
}
sub c {
my $x = shift;
my $result = 0;
for my $k (0..$maxn-1) {
$result += (-1)**$k * &p($x, 2 * $k + 0) / &factorial(2 * $k + 0);
}
return $result;
}
sub s {
my $x = shift;
my $result = 0;
for my $k (0..$maxn-1) {
$result += (-1)**$k * &p($x, 2 * $k + 1) / &factorial(2 * $k + 1);
}
return $result;
}
# ($x - 0)($x - 1)(...)($x - ($n-1))
sub p {
my ($x, $n) = @_;
my $result = 1;
for my $k (0..$n-1) {
$result *= $x - $k;
}
return $result;
}
sub factorial {
my $n = shift;
my $result = 1;
while ($n > 0) {
$result *= $n;
$n--;
}
return $result;
}
\end{verbatim}
動かしてみたところ、次のようになりました。
\begin{verbatim}
0 1 0
1 1 1
2 0 2
3 -2 2
4 -4 0
5 -4 -4
6 0 -8
7 8 -8
8 16 0
9 16 16
10 0 32
11 -32 32
12 -64 0
13 -64 -64
14 0 -128
15 128 -128
\end{verbatim}
何と。
「離散系バージョンの関数探し」
の$f(x)$および$g(x)$と同じ結果ではありませんか。
\par
同じ差分方程式でおなじ初期値なのだから、当然なのか…。
くすん。
\newpage
\begin{thebibliography}{9}
\bibitem{diffsum}
結城浩,
「ミルカさんの隣で」,
\verb|http://www.hyuki.com/story/diffsum.html|,
2005.
\bibitem{diffsum2}
結城浩,
「離散系バージョンの関数探し」,
\verb|http://www.hyuki.com/story/diffsum2.html|,
2005.
\bibitem{present}
吉田武,
『オイラーの贈り物』(p.~216),
ちくま学芸文庫,
ISBN4-480-08675-7,
2001.
\end{thebibliography}
\end{document}
追記
(み)さんから
「三角関数の離散系バージョン」の話ですが,
s0(x) = 0, 1, 1, 0,-1,-1, 0, ... (x=0,1,2,3,4,5,6,...) s1(x) = 1, 0,-1,-1, 0, 1, 1, ... s2(x) =-1,-1, 0, 1, 1, 0,-1, ...
というのをたまたま思いつきました. D sin(x) = sin (x + pi/2) に対し,この s0(x) は Delta s0(x) = s0(x + 1) という性質があるので,s0(x) は sin x の離散系バージョンの一つといってよいかと思います.
結城から
なるほど。cos, sinのデュオではなく、 s0, s1, s2のトリオなのですね。 面白いです!
ぴろさんから
f(0) = 1 g(0) = 0 f(x+1) = (f(x) + g(x)) / sqrt(2) g(x+1) = (g(x) - f(x)) / sqrt(2)
でどうでしょう?
結城から
なるほど。 確かにこれだと、周期性をもった関数になりますね。 ありがとうございます。 ただ…、これだと「三角関数の微分方程式に対応する差分方程式から導いたよ」 という感覚が薄れてしまうのでちょっと残念です。
ぴろさんから、さらに
f(x+2) = f(x+1) + g(x+1) = 2g(x)
g(x+2) = g(x+1) - f(x+1) = -2f(x)
から三角関数の加法定理を使って
f(x+1)g(1)+g(x+1)f(1) = 2g(x)
f(x+1)f(1)-g(x+1)g(1) = 2f(x)
ここでf(1)=g(1)とするとf(1)^2+g(1)^2=1からf(1)=g(1)=sqrt(2)で
f(x+1)=(f(x)+g(x))/sqrt(2)
g(x+1)=(g(x)-f(x))/sqrt(2)
ということで差分方程式から導いたっぽくないでしょうか?
結城から
うんうん、導いたっぽいですね! (^_^)V
以下、そのほかの方々からのリンクです。
さらにリンクを追記(2005-05-14)
「ミルカさんの隣で」を書いてから、 LaTeXでEulerフォントで数式を出すのが楽しくなり、 その勢いで今度は 「離散系バージョンの関数探し」という文章を書きました。
Eulerフォント( AMS Euler) というのは、Hermann Zapf(ザフ)がデザインした数式用のフォントです。 『コンピュータの数学』で使われており、同書の著者前書きには、 「ザフの設計の基本方針は, 字の上手な数学者が手書きした場合の数学の薫りをもたせることである.」(p.viii) と書かれています。
確かにEulerフォントは、
全体的に手書きの風合いがあります(たとえば数字のゼロの上がとんがっている)。
このフォントで数式を出力すると、自分が書いたものでも、
まるで偉大な数学者の手によるもののような気分になれます(^_^)。
Eulerフォントを使うときのLaTeXの書き方は、たとえばこんな感じです。 途中に \usepackage{euler} と書いてあるのがわかると思います。
\documentclass{jsarticle}
\title{ミルカさんの隣で\footnote{http://www.hyuki.com/story/diffsum.html}}
\author{結城浩\footnote{http://www.hyuki.com/\ Hiroshi Yuki \copyright\ 2005, All rights reserved.}}
\date{2005年5月}
\usepackage{amsmath}
\usepackage{euler}
\begin{document}
\maketitle
...
\end{document}
さて、 「ミルカさんの隣で」は連続と離散の話(微分と差分の話)でした。 あの文章にはx^nしか出てきませんでしたが、 ほかの関数についても少し調べてPDFに書いてみました。 ご興味のある方はお読みください。 また、誤りのご指摘、ご意見やご感想もぜひお送りください。
追記
「Hermann Zapfはヘルマン・ツァップと読むはず」というご指摘をいただきましたが、 『コンピュータの数学』では「Hermann Zapf(ザフ)」と表記されていましたので、 それに倣っています。
今日は母の日ですね。 夕食後、実家の母に電話をしておしゃべりをしました。 おかあさん、ありがとう。
『プログラマの数学』の読者さんからうれしいお便りをいただいたので、許可を得てご紹介します。
読者さんから
GW中を利用して「プログラマの数学」を読んでいます。 まるっきりの文系人間で、中学生くらいの頃から数学が苦手だったのですが、 この本のおかげで「もしかして数学が嫌いだったわけじゃなくて、学校の授業に ついていけなかったから嫌いと思っていただけで、ホントは好きなんじゃないか」 と思うほど、面白く読ませてもらっています。
結城から
ご愛読ありがとうございます。 また、とても励みになるフィードバックも感謝です。 拙著を楽しんで読んでいただけているようで、とてもうれしいです!
物語「ミルカさんの隣で」を公開します。 ぜひ、ご感想をお寄せください。
まずは 薬品格納クイズ(問題編)をお読みください。
結城の解答
縦に、
A B
と並んでいる状態を、 (AB) と表記することにします。またたとえば、
ABA BAA
のように並んでいる状態を、 (AB-BA-AA) と表記することにします。
まずは、列数nが1, 2の場合で数えます。
n=1のときは、以下の3通り。
(AA) (AB) (BA)
n=2のときは、以下の7通り。
(AA-AA) (AA-AB) (AA-BA) (AB-AA) (AB-BA) (BA-AA) (BA-AB)
n=3のときは、以下の17通り。
(AA-AA-AA) (AA-AA-AB) (AA-AA-BA) (AA-AB-AA) (AA-AB-BA) (AA-BA-AA) (AA-BA-AB) (AB-AA-AA) (AB-AA-AB) (AB-AA-BA) (AB-BA-AA) (AB-BA-AB) (BA-AA-AA) (BA-AA-AB) (BA-AA-BA) (BA-AB-AA) (BA-AB-BA)
ここまで作ってくると、次のことがわかります。
n列で、右端がAAになる場合の数をa(n), 右端がABになる場合の数をb(n), 右端がBAになる場合の数をc(n)と置くことにします。 すると、上で調べたことから、次のことがわかります(n >= 1)。
a(1) = 1 b(1) = 1 c(1) = 1 a(n+1) = a(n) + b(n) + c(n) b(n+1) = a(n) + c(n) c(n+1) = a(n) + b(n)
対称性からb(n)=c(n)が成り立ちます。式を整理すると、
a(n+1) = a(n) + 2b(n) b(n+1) = a(n) + b(n)
求める値はa(6) + b(6) + c(6)なので、あとは順番に計算します。
a(6) + b(6) + c(6) = a(6) + 2b(6)
= (a(5) + 2b(5)) + 2(a(5) + b(5))
= 3a(5) + 4b(5)
= 3(a(4) + 2b(4)) + 4(a(4) + b(4))
= 7a(4) + 10b(4)
= 7(a(3) + 2b(3)) + 10(a(3) + b(3))
= 17a(3) + 24b(3)
= 17(a(2) + 2b(2)) + 24(a(2) + b(2))
= 41a(2) + 58b(2)
= 41(a(1) + 2b(1)) + 58(a(1) + b(1))
= 99a(1) + 140b(1)
= 99 + 140
= 239
答え:239通り
以下に、読者さんからの解答から何通か紹介します(すべてではありません)。 順番はほぼ到着順です。 みなさん、ご解答ありがとうございます。
ぞりおさんの答え
はじめまして。 クイズやってみました。
薬品箱の状態を12ビットで表します。 あらかじめ爆発するケースを配列で持っておきます。 0から0xFFFFFFまでの数を、爆発するケースと比較して 判定しました。 答えは239通りでした。
西尾泰和さんの答え
notogawaさんの答え
上下にAA,AB,BAと並ぶパターンをそれぞれα、β、γ、とする
βとγは連続して並ぶことの無い並べ方がここでは適当である
(例えば左から)順番にα、β、γを適当に並べていって
n番目にα、β、γが配置される並べ方の数をそれぞれ
α_n、β_n、γ_nとするとα_1=β_1=γ_1=1で、
α_{n+1}=α_n+β_n+γ_n (αはどの後でもOK)
β_{n+1}=α_n +γ_n (ββは不適当)
γ_{n+1}=α_n+β_n (γγは不適当)
であり、このとき欲しい解はα_6+β_6+γ_6=α_7である
簡単のためS_n=α_n、T_n=β_n+γ_nとするとS_1=1、T_1=2で
S_{n+1}=S_n+T_n
T_{n+1}=2S_n+T_n
となりこれらよりSについての三項間漸化式は
S_{n+2}=2S_{n+1}+S_n (S_1=1, S_2=S_1+T_1=3)
この後S_nの一般項を求めてもいいが
ルートとか出てきてかえって面倒で、
欲しい値はα_7=S_7なので順番に
S_3=2S_2+S_1=7
S_4=2S_3+S_2=17
S_5=2S_4+S_3=41
S_6=2S_5+S_4=99
S_7=2S_6+S_5=239
以上より239通り
ishigakiさんの解答
=comment
こんばんは。何度かPerlクイズなどでメールを差し上げた
石垣です。ぐずぐずしているうちに答えが出てしまいまし
たが、上下左右の反転チェックまで考慮したものはないよ
うでしたので送ってみます。ちなみに上下左右固定したも
のは239通り、反転したものを同一とみなすと66通りでした。
以下、チェックに利用したスクリプトのソースです。
=cut
#!/usr/bin/perl
# ここでは隣接すると発火してしまうBをどのように
# 置くかだけが問題なので、まずは横の列に注目して、
# Bが隣接しない組み合わせを抽出します。ビット演
# 算で1マス分シフトしたとき、ビット積が0になれば
# 隣接したBはない、と判断できます。
my $col = 6; # 横の列数
my $max = 2 ** $col;
my @row = ();
for(my $ct = 0; $ct < $max; $ct++) {
my $shifted = $ct >> 1;
next if $shifted & $ct;
# ビット演算をクリアした組み合わせを保存
push @row, $ct;
}
# 今度は縦の列について考えます。横の列の場合と同
# 様に、上下のビット積が0になれば隣接したBはない
# と判断できます。
my $total = 0;
my %exists = ();
foreach my $top (@row) {
foreach my $bottom (@row) {
# 上下のビット積が0でなければ無視します
next if $top & $bottom;
# 以下の五行は参考までに。
# コメントを外すと、上下左右を反転したと
# きに同じ形になる組み合わせを排除します
# next if $exists{ bin($top).bin($bottom) };
# $exists{ bin($top).bin($bottom) } = 1;
# $exists{ bin($bottom).bin($top) } = 1;
# $exists{ rev($top).rev($bottom) } = 1;
# $exists{ rev($bottom).rev($top) } = 1;
$total++;
}
}
print "total: $total\n";
# 以下、サブルーチン
# bin は $col 桁数の2進数を、rev の方は bin で得ら
# れた2進数の左右反転したものを返します。
sub bin { substr(unpack('B16',pack('n',shift)),-1 * $col); }
sub rev { join('',reverse(split('',sprintf('%s',bin(shift))))); }
Knuth先生によるThe Art Of Computer Programming, Volume 4はまだ出ていませんけれど、 いくつかの章はすでにFascicle(ファスィクル, 分冊)という名の下に出版されていますのでご紹介します。
ちなみに、その続きの章はPre-Fascicleとして、 Knuth先生のページからダウンロードできます。
ついでに、日本語版が刊行されているVolume 1, 2もリンクしておきます。 たしかVolume 3の日本語版も2005年には刊行されるはず。
antを使って、${builddir}にあるクラスファイルと、 ${resourcedir}にあるJPEGファイルを ひとつのjarにまとめようと思った。 最初は、次のようにした。
<jar jarfile="${jardir}/${ant.project.name}-${DSTAMP}.jar" basedir="${builddir}" />
<jar jarfile="${jardir}/${ant.project.name}-${DSTAMP}.jar" basedir="${resourcedir}" />
でも、jar tvfで見てみると、JPEGファイルしか入っていない。あれれ? そうか、2回目のjarタスクで、新しく作り直しちゃったんだね。
antのマニュアルを見直して、次のように修正。
<jar destfile="${jardir}/${ant.project.name}-${DSTAMP}.jar">
<fileset dir="${builddir}" />
<fileset dir="${resourcedir}" />
</jar>
CodeZine(コードジン)という「開発者のための実装系Webマガジン」に、 「モジュールを使わないシンプルなアクセスカウンタ」という記事を書いてみました。
まあ要するに、Perlで作ったアクセスカウンタなのですが、 GDなどのモジュールを使わず、 外部のGIFファイルも使わないというのが特徴です。 ご興味のある方はお読みください。
以下のように6×2=12個の区画に分かれている薬品箱があります。 ここにA, Bという二種類の薬品を格納します。 ところが、薬品Bは縦や横に隣り合った区画に並べて格納すると発火してしまいます(斜めなら大丈夫。薬品Aは並べても大丈夫)。 このとき、発火しない格納方法は全部で何通りあるでしょうか。
薬品箱
□□□□□□ □□□□□□
発火する格納例
ABABAB AABAAB ABABBB AABAAA
発火しない格納例
ABABAB AABAAA AAAABA BAAAAA AABAAA BAAABA
解答したい方はfeedbackからどうぞ。 解答は、 日記などで公開させていただくかもしれませんので、 そのおつもりで。 もちろん、自分のblogに書いて、URLを送ってくださってもOKです。
Enjoy!
なお、この問題は 『マスター・オブ・場合の数』(p.24)によります。 これはとても楽しい本で、拙著 『プログラマの数学』の中でも何点か参考にさせてもらっています。