続々・マーデルング定数（ＣｓＣｌの場合）

　ＣｓＣｌ型構造では，１個のＣｓイオンに最も接近しているのは８個のＣｌイオンである．この最近接距離をｒとすると，第２近接は（２／√３）ｒの距離にある６個のＣｓイオン，第３近接は（２√２／√３）ｒの距離にある１２個のＣｓイオン，第４近接は（√１１／√３）ｒの距離にある２４個のＣｌイオン，さらに８個のＣｓイオンが２ｒの距離（第５近接）に存在する．

　したがって，ＣｓＣｌのマーデルング級数は，

　　８－６／（２／√３）－１２／（２√２／√３）＋２４／（√１１／√３）－８／４－・・・

になる．同様に，ＮａＣｌのマーデルング級数は，

　　６－１２／√２＋８／√３－６／２＋２４／√５－・・・

であり，ＮａＣｌとＣｓＣｌではイオンの並び方が異なっている．

　ＣｓＣｌのマーデルング級数を，実際に計算すると

　Ｍ＝８－５．１９６－７．３４８＋１２．５３４－４－・・・

しかし，この方法では収束はよくない．マーデルング級数の困難な点は，級数のなかの陽イオンの項と陰イオンの項が別々に発散することである．すなわち，∞－∞の不定形であり，条件つき収束のため，マーデルング級数の収束はきわめて遅く，プログラム的には収束しないという結論になってしまう．→【補】

　上記の求め方はいわば極座標的な級数和であって，座標軸を変えることによって，マーデルング級数の足し上げの順序が変わってくる．マーデルング級数の面白い点は，順番を替えて足し合わせると収束させることができることである．しかも，プログラムの組み方によってはノロノロ這い回るカメをピョンピョン跳ね回るウサギに転化させることさえ可能なのである．

　実のところ，この一連のコラムを書くまでは，順番を入れ替えると収束する級数というのは数学上の病理学的な級数であって，応用分野では存在しないであろうと思っていた．小生にとってマーデルング級数の出現は意外な，青天の霹靂というべき事態であった．

　さて，前回掲げたプログラムで，ＣｓＣｌのマーデルング定数が正しい値に収束しなかった原因と対策について，東北大学・工学部・材料化学の遠藤忠教授よりご教示いただいたところ，ＮａＣｌではアニオンとカチオンが交互に並び，格子面が電気的に中性になるが，ＣｓＣｌでは電気的中性条件を満たさないことが原因とのことであった．　その対応策は如何に・・・・・？

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

　面心立方格子や対心立方格子では，直交座標を基本としている．直交座標軸は空間中の点の位置を表すのに最も取り扱いが簡単である．しかし，３つベクトルａ↑，ｂ↑，ｃ↑の選び方は一義的には決まらず，いろいろな選び方がある．そこで，３つのベクトルをそれぞれの構造の単純並進ベクトルと呼ばれるものに採ってみた．そうすると，対心立方格子ではａ＝ｂ＝ｃ，α＝β＝γ＝１０９．４９°，面心立方格子ではａ＝ｂ＝ｃ，α＝β＝γ＝６０°の菱形体格子に変換されることになる．

　これによって，ＣｓＣｌでは，格子点に交互にカチオンとアニオンが並ぶ菱形体格子が基本単位になり，これを各軸方向に２倍して８個の集まりを考えると，ＮａＣｌ型の単位格子と同様のイオン配置をした大きな菱形体格子ができあがる．非直交座標軸に変換されるのだが，この菱形体格子ではａ＝ｂ＝ｃ，α＝β＝γ＝１０９．４９°であることに注意して計算してみると，実際にマーデルング定数は収束値（１．７６２６７）に向かった．以下に，そのプログラムと計算結果を示す．

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

（プログラム）

1010 ' *** Madelung constant ***

1040 PRINT "dimension=";DIMENSION

1050 INPUT "distance =";D

1110 PRINT "Madelung=";SUM

1370　 FOR X=-D TO D-1

1380　　 FOR Y=-D TO D-1

1390　　　 FOR Z=-D TO D-1

1400　　　　 PRINT "*";

1410　　　　 GOSUB *VERTEX:UV=U/8

1420　　　　 GOSUB *EDGE:UE=U/4

1430　　　　 GOSUB *FACE:UF=U/2

1440　　　　 GOSUB *CELL:UC=U

1450　　　　 SUM=SUM-UV+UE-UF+UC

1460　　　 NEXT Z

1480　　 PRINT X,SUM

1540 S1=XP*XP*AA*AA+YP*YP*BB*BB+ZP*ZP*CC*CC

1550 S2=2*(XP*YP*AA*BB*AB+YP*ZP*BB*CC*BC+ZP*XP*CC*AA*CA)

1560 R2=ABS(S1+S2):R=SQR(R2)

1570 IF R<=0 THEN 1610

1680 XP=X+1:YP=Y:ZP=Z

1710 XP=X:YP=Y+1:ZP=Z

1740 XP=X+1:YP=Y+1:ZP=Z

1770 XP=X:YP=Y:ZP=Z+1

1800 XP=X+1:YP=Y:ZP=Z+1

1830 XP=X:YP=Y+1:ZP=Z+1

1860 XP=X+1:YP=Y+1:ZP=Z+1

1920 XP=X+.5:YP=Y:ZP=Z

1950 XP=X:YP=Y+.5:ZP=Z

1980 XP=X+1:YP=Y+.5:ZP=Z

2010 XP=X+.5:YP=Y+1:ZP=Z

2040 XP=X:YP=Y:ZP=Z+.5

2070 XP=X+1:YP=Y:ZP=Z+.5

2100 XP=X:YP=Y+1:ZP=Z+.5

2130 XP=X+1:YP=Y+1:ZP=Z+.5

2160 XP=X+.5:YP=Y:ZP=Z+1

2190 XP=X:YP=Y+.5:ZP=Z+1

2220 XP=X+1:YP=Y+.5:ZP=Z+1

2250 XP=X+.5:YP=Y+1:ZP=Z+1

2310 XP=X+.5:YP=Y+.5:ZP=Z

2340 XP=X+.5:YP=Y:ZP=Z+.5

2370 XP=X:YP=Y+.5:ZP=Z+.5

2400 XP=X+1:YP=Y+.5:ZP=Z+.5

2430 XP=X+.5:YP=Y+1:ZP=Z+.5

2460 XP=X+.5:YP=Y+.5:ZP=Z+1

2520 XP=X+.5:YP=Y+.5:ZP=Z+.5

　前回のプログラムと異なっている点は，

１）サブルーチン＊ＣｓＣｌで，塩化セシウムの格子定数を与えている．変数ＡＡ，ＢＢ，ＣＣはそれぞれ格子定数ａ，ｂ，ｃに，変数ＡＢ，ＢＣ，ＣＡは格子定数γ，α，βに対応している．塩化セシウムの場合，等長・等角の格子に設定してあるが，もちろん，これらの制限を緩めることも可能で，任意の平行六面体であってよい．

２）非直交座標における原点からの距離の計算（１５４０行～１５６０行）

３）直交座標ではあれば，対称性が高いので，第１象限だけを数え上げてその結果を全空間に外挿することができる．その場合は約８倍の時間節約効果が期待できるが，任意の座標軸を取り扱う場合，２）の距離計算は原点に対する対称性しか成り立たたないから，その節約効果は高々２倍にすぎない．そのため，サブルーチン＊ＳＵＭ３の因子ＭをＭ＝１とし特別な加速操作を行っていないこと（１５９０行）

だけであり，プログラムの書き換えもたった数行で済んでいる．

　また，１６３０行～２５４０行で原点の位置を対心から頂点に変更しているが，中性（もしくはそれに近い）領域を１つの細胞とみなし，原点から外に向かって順々に各細胞ごとに加えてゆくやり方や細胞境界に存在するイオンに対しては分割して分数の電荷を与えるやり方など，基本的な骨格に変更はない．

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

（結果）単精度で計算した結果を示す．なお，加算される細胞数は

　　８個→６４個→２１６個→８ｄ^3→・・・となる．

ａ）ＣｓＣｌ結晶の場合

　　ｄ＝１　→　α＝１．７９１０　（細胞数：８）

　　ｄ＝２　→　α＝１．７６３０　（細胞数：６４）

　　ｄ＝３　→　α＝１．７６２７　（細胞数：２１６）

　　ｄ＝４　→　α＝１．７６２７　（細胞数：５１２）

　　ｄ＝５　→　α＝１．７６２７　（細胞数：１０００）

　　ｄ＝６　→　α＝１．７６２７　（細胞数：１７２８）

　　ｄ＝７　→　α＝１．７６２７　（細胞数：２７４４）

　　ｄ＝８　→　α＝１．７６２８　（細胞数：４０９６）

　このように，非直交座標に変換して空間を切り直すとあっという間に正しい値：α＝１．７６２７６に収束した．ＣｓＣｌの場合，細胞全体としては中性であっても，格子面が中性条件を満たさないことは小生も気づいていたのであるが，直交座標を用いる限り，この点は解消されないし正しい値も求まらない．

　直交ベクトルはわかりやすいし，便利なものであるから慣用的に用いられているのであろうが，ＣｓＣｌの場合は確かに直交座標があだとなってしまう．ベクトルの選び方・格子変換が収束に対してこれほど有用だとは思いもよらず，驚異的でさえあった．

ｂ）ＮａＣｌ結晶の場合（菱形面格子）

　１２００行～１２６０行で格子定数を設定

　　ｄ＝１　→　α＝１．７１５９　（細胞数：８）

　　ｄ＝２　→　α＝１．７４７８　（細胞数：６４）

　　ｄ＝３　→　α＝１．７４７６　（細胞数：２１６）

　　ｄ＝４　→　α＝１．７４７６　（細胞数：５１２）

　　ｄ＝５　→　α＝１．７４７６　（細胞数：１０００）

　　ｄ＝６　→　α＝１．７４７６　（細胞数：１７２８）

　　ｄ＝７　→　α＝１．７４７６　（細胞数：２７４４）

　　ｄ＝８　→　α＝１．７４７６　（細胞数：４０９６）

ｃ）ＮａＣｌ結晶の場合（面心立方格子）

　１１４０行～１１８０行で格子定数を設定

　　ｄ＝１　→　α＝１．７５１８　（細胞数：８）

　　ｄ＝２　→　α＝１．７４７７　（細胞数：６４）

　　ｄ＝３　→　α＝１．７４７６　（細胞数：２１６）

　　ｄ＝４　→　α＝１．７４７６　（細胞数：５１２）

　　ｄ＝５　→　α＝１．７４７６　（細胞数：１０００）

　　ｄ＝６　→　α＝１．７４７６　（細胞数：１７２８）

　　ｄ＝７　→　α＝１．７４７６　（細胞数：２７４４）

　　ｄ＝８　→　α＝１．７４７６　（細胞数：４０９６）

　また，ＮａＣｌについてはｃ）のように，直交座標上でもマーデルング定数が求まっていたのだが，ｂ）のように互いに６０°で交わる菱面体格子（α＝β＝γ＝６０°）で表してもα＝１．７４７５６に収束し，直交軸・非直交軸のいずれからも同じ結果が得られたことになる．

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

【補】絶対収束・条件収束

　Σａnの絶対値級数Σ｜ａn｜が収束するとき，Σａnは絶対収束するという．また，Σａnは収束するが，Σ｜ａn｜は収束しない級数は，条件収束するという．

　マーデルングの交代級数では，１種のディリクレ級数→【補】

　　Σ±ａn／ｎ^s　　（ｓ＝１／２）

の計算をしていることになる．数学に堪能な方であれば，この級数が絶対収束するかどうか調べるべきところであるが，小生，巧い判定法を思いつかない．しかし，結晶の周期的構造の物理的意味を考えると絶対収束しないことは明らかであろう．

　　１／１－１／２＋１／３－１／４＋・・・

　　１／１＋１／２＋１／３＋１／４＋・・・→　＋∞

の交代級数である．この値は対数関数のマクローリン展開

　　ｌｏｇ（１＋ｘ）＝ｘ－１／２ｘ2 ＋１／３ｘ3 －１／４ｘ4 ＋・・・

においてｘ＝１とおくとｌｏｇ２に収束することがわかる．この値はメルカトールの定数とかグレゴリーの定数と呼ばれている定数である．

　ところで，交代級数では，元の級数の項の順番を変えると収束値が変動してしまう．たとえば，負項を正項に変えて，あとでその２倍を引くと，

　１／１－１／２＋１／３－１／４＋・・・

＝（１／１＋１／２＋１／３＋１／４＋・・・）－２（１／２＋１／４＋１／６＋１／８＋・・・）

＝（１／１＋１／２＋１／３＋１／４＋・・・）－（１／１＋１／２＋１／３＋１／４＋・・・）

　また，この交代級数は奇数の逆数と偶数の逆数に－１をかけたものからできているが，足し合わせる順序が違う級数，たとえば，負の項が２つの連続する正の項をはさんで現れる級数：

　　｛１／１＋１／３－１／２｝＋｛１／５＋１／７－１／４｝＋・・・

では３／２ｌｏｇ２に収束する．また，正の項に引き続いて負の項が２つの連続する級数：

　　｛１／１－１／２－１／４｝＋｛１／３－１／６－１／８｝＋・・・

は１／２ｌｏｇ２に収束することがわかっている．

　　｛１／１－１／２－１／４｝＋｛１／３－１／６－１／８｝＋・・・

　　＝１／２ｌｏｇ２を示す．

　与えられた級数は

　Σ｛１／（２ｎ－１）－１／２（２ｎ－１）－１／（２（２ｎ－１）＋２）｝

＝Σ｛１／（４ｎ－２）－１／４ｎ｝

　一方，１／１－１／２＋１／３－１／４＋・・・＝ｌｏｇ２より

１／２ｌｏｇ２＝１／２－１／４＋１／６－１／８＋・・・

　　　　　　　＝（１／２－１／４）＋（１／６－１／８）＋・・・

　　　　　　　＝Σ｛１／（４ｎ－２）－１／４ｎ｝

　これらは無限のパラドックスの一つの例である．有限級数ならば，足し算の順序に入れ替えは自由にできるが，無限級数となると話はまったく違ってくる．正の項と負の項がいずれも絶対収束するとき，級数の和の順番は勝手に変えてもよいのであるが，そうでない場合は，足す順序によっては級数の和が異なってくる．実は，条件収束級数の場合，級数の項の順番を適当に変えるとどんな値にでも収束させることができることが知られている．

定理（１）：絶対収束級数は項の順序をどのように変えても絶対収束し，和も変わらない．（ディリクレ）

定理（２）：条件収束級数は項の順序を適当に変えれば，指定された値（±∞でもよい）を和にもつようにも，振動するようにもできる．（リーマン）

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

【補】ディリクレ級数とメリン変換

　与えられた数列｛ａn｝の性質を知ろうとしたとき，数列｛ａn｝を直接調べるのではなく，たとえば，ｘを変数とする級数関数

　　φ（ｘ）＝ａ0 ＋ａ1 ｘ＋ａ2 ｘ2 ＋・・・＋ａn ｘn＋・・・

をつくり，この関数を観察します．数列の性質は何らかの形で関数の振る舞いに反映していると考えられます．こうして数列｛ａn ｝の数論的性質（離散的）を，関数φ（ｘ）の解析的性質（連続的）によって解明していこうとする発想が母関数のアイディアです．

　　φ（ｘ）＝Σａnｋ（ｎ，ｘ）

において，関数ｋ（ｎ，ｘ）をカーネル（核関数）と呼びます．通常型母関数はカーネルとしてk(n,x)=x^nを，指数型母関数ではk(n,x)=x^n/n!を利用しています．k(n,x)=1/n^xを利用した場合，その母関数をディリクレ母関数と呼びます．｛１，１，１，・・・｝のディリクレ母関数はゼータ関数ζ(x)=Σ1/n^xです．

　また，関数f(x)に対して，積分

　　h(s)=integral(0-∞）x^(s-1)f(x)dx

が存在するとき，これを関数f(x)のメリン変換といいます．exp(-x)のメリン変換はガンマ関数Γ(s)であることより，ガンマ関数の定義も１種のメリン変換です．

　　Γ(s)=integral(0-∞）x^(s-1)exp(-x)dx

　ガンマ関数の定義式より

　　integral(0-∞）x^(s-1)exp(-nx)dx=Γ(s)n^(-s)

ですから，ディリクレ級数Σann^-sについて

　　Σann^-s=1/Γ(s)integral(0,∞)(Σanexp(-nt))t^(s-1)dt

が得られます．この式はディリクレ級数f(s)=Σann-sと同じ係数をもつベキ級数F(z)=Σanz^nはメリン変換

　　f(s)=1/Γ(s)integral(0,∞)F(exp(-t))t^(s-1)dt

によって互いに結ばれていることを意味します．

ζ(ｓ)=Σ1/ｎ^sにおいて，F(exp(-t))=Σexp(-nt)=1/(exp(t)-1)

φ（ｓ）＝Σ(-1)^(n-1)/n^s=(1-2^(1-s))ζ(ｓ)において，

F(exp(-t))=Σ(-1)^(n-1)exp(-nt)=1/(exp(t)+1)

Ｌ（ｓ）＝　１／１s －１／３s ＋１／５s －１／７s ＋・・・

において，F(exp(-t))=1/(exp(t)+exp(-t))

Γ(ｓ)ζ(ｓ)=integral(0-∞)x^(s-1)/(e^x-1)dx

Γ(ｓ)ζ(ｓ)(1-2^(1-x))=integral(0-∞)x^(s-1)/(e^x+1)dx

Ｌ（ｓ）=1/Γ(s)integral(0,∞)t^(s-1)/(exp(t)+exp(-t))dt

　関数論において，ベキ級数は基本的な役割を演じますが，メリン変換によって，ベキ級数の性質からディリクレ級数の性質を導いたり，その逆も可能になります．ゼータ関数は最も簡単かつ最も重要なディリクレ級数f(s)=Σann-sですが，メリン変換は（解析的数論における）ゼータ関数と（関数論における）保型関数（ある種の２重周期的挙動をする複素変数関数）を結ぶ装置として，数論の世界では決定的に重要な意味をもっています．

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝