固有値計算の並列化とマクロタスク機能

奈良県立医科大学　平井　國友

　スーパーコンピュータＳＸ－４での固有値計算の並列化とマクロタスク機能について報告する。これまでの報告では、主に自動並列化について調べ、固有値計算の並列化がＳＸ－４で有効に機能していることを結論として述べた。ＳＸ－４のＦＯＲＴＲＡＮコンパイラの自動並列化はかなり有効で、ワークステーションなどで使っている計算プログラムコードをほとんど変更しなくても十分に並列化の恩恵を得られる。これはＳＸ－４の非常に優れた点といえる。ただし、固有値計算の場合には、ＡＳＬライブラリ中に実対称行列の全固有値・全固有ベクトルを求めるＱＣＳＭＡＡ（ＤＣＳＭＡＡの並列化版）サブルーチンが用意されていることが並列化にかなり有利に働いている。このため、並列版（マニュアル「ＡＳＬ/ＳＸ利用の手引＜並列処理機能編＞」を参照）が用意されていない場合には、並列化の恩恵がそれほど大きくない可能性もある。

　さて、これまでは、行列1個の固有値計算について論議してきたといえる。しかし、現実問題としては、多数の行列の固有値をまとめて計算することのほうが多いと思われる。このような場合、行列1個の固有値計算を自動並列化すなわちマイクロタスク機能を使って並列化するよりは、異なる行列（行列の大きさは同じ）の固有値計算を並列化して行なうほうが効率が良い可能性が高い。そこで、今年度はこのような並列化、すなわち、マクロタスク機能を使った並列化について報告する。マクロタスク機能については、センターニュースではほとんど説明されていないようなので、不案内なユーザも多いであろう。そこで簡単な使用例として紹介する。前年度は解析ツールの一つであるＡＮＡＬＹＺＥＲの使用法についても紹介したが、そのときと同様に、必要な情報がマニュアルのどこに記述されているかを探し当てるのがそれほど容易ではなく、マクロタスク機能を何とか使えるまでにかなりの時間を費やした。ここで紹介する使用例が、多くのユーザーの参考になれば幸いである。

　最初に、マイクロタスク機能とマクロタスク機能の違いについて簡単に説明する（詳細についてはマニュアル「ＦＯＲＴＲＡＮ77/ＳＸ並列処理機能利用の手引」を参照）。両者は共に並列化のための機能ではあるが、マイクロタスク機能が”ＤＯループ単位の並列化”であるのに対して、マクロタスク機能は”サブルーチン単位の並列化”であるという違いがある。マイクロタスク機能は自動並列化あるいは各自が並列化指示行を追加することで利用でき比較的容易に並列化の恩恵を得ることができる。しかし、問題なのはマイクロタスク機能ではＤＯループ中にサブルーチンがあると並列化できないという点である。このためマイクロタスク機能ではプログラムの局所的な部分でしか並列化できないということになる。一方、マクロタスク機能では特別な組み込みサブルーチンを使って利用できる（したがって、プログラムコードをある程度変更する必要がある）。マクロタスク機能は大きな単位の並列化が可能であり、一般的な意味での並列化の機能といえる。ところが、上で述べたように、このマクロタスク機能に関する説明は十分であるとは言い難い。この点では今後改善すべきであろう。

　マクロタスク機能を利用するためのプログラムコードの変更は、通常、２段階（あるいは３段階）で行なうことになる。固有値計算を例にすると、第１段階では、行列要素の計算とその固有値計算の部分を１個のサブルーチン（ここでは、ＥＩＧＥＮという名前にする）にまとめてしまう。すなわち、ｃａｌｌ　eigen（ｉ）とすると、ｉ番目（　i＝１，２，．．．，N）の行列の固有値計算ができるようにする（Nは計算する行列の個数）。このときサブルーチンの引数が番号を指定する i　だけになるように、COMMON文などを使って変更する。サブルーチンの引数が他にあっても差し支えないが、第１段階としては、まず、この形に変更するのが確実な方法であると考えられる。
この変更で問題がなければ、第２段階では、以下のようなに、組込みサブルーチンを使用して並列化を行う。

　ｐｒｏｇｒａｍ　ｍａｉｎ
　ｅｘｔｅｒｎａｌ　ｅｉｇｅｎ
　ｇｌｏｂａｌ　ｃｏｍｍｏｎ/ｔａｓｋ/ｉｔｖａｒ（４），ｉｔｐｒａｍ（４），ｎｏ（４）
　　　．．．．
　ｃａｌｌ　ｐｓｔｕｎｅ（’ｍａｘｃｐｕ’，４）
　ｄｏ　２　ｉ＝１，４
２ｎｏ（ｉ）＝ｉ
　ｄｏ　２　ｉ＝２，４
４ｃａｌｌ　ｐｔｆｏｒｋ（ｉｔｖａｒ（ｉ），ｉｔｐｒａｍ（ｉ），ｅｉｇｅｎ，ｎｏ（ｉ））
　ｃａｌｌ　ｅｉｇｅｎ（ｎｏ（１））
　ｃａｌｌ　ｐｔｊｏｉｎ（ｉｔｖａｒ（２），ｉｔｖａｒ（３），ｉｔｖａｒ（４））
　　　．．．．
　ｓｔｏｐ
　ｅｎｄ

　ここでは、Ｎ＝４の場合を考え、４個の行列の固有値計算を組込みサブルーチンＰＳＴＵＮＥ、ＰＴＦＯＲＫ、ＰＴＪＯＩＮを用いて並列化している。サブルーチンＰＳＴＵＮＥは並列化に関するパラメータを変更するもので、ここでは、ＣＰＵ数を４にしている（その他のパラメータについてはマニュアルを参照）。サブルーチンＰＴＦＯＲＫで並列実行を開始し、 i ＝２，３，４番目の行列の固有値計算をもとのＣＰＵとは別の３つのＣＰＵに振り分けて実行させ、もとのＣＰＵでは i ＝１番目の行列の固有値計算を実行する（call　eigen(no(1))）。そしてサブルーチンPTJOINで各CPUの計算が終了するのを待ち合わせる。ここで、サブルーチンPTFORKの４番目以下の引数はサブルーチンEIGENの引数である。整数型の変数itvar，itpramは必須であり。通常この例のようにしておく。整数型の変数noは３つのCPUに i ＝２，３，４の値を同時に振り分けるために必要であり no( i )= i の文も含めて、実質上、これも必須といえる。これらの変数はGLOBAL　COMMON文で共通データに割り付けておく。また、サブルーチン名EIGENがサブルーチンPTFORKの３番目の引数になっているので、EXTERNAL文で宣言しておかなければならない。
　並列化のための変更はとりあえずこれで十分で、後は、ｆ７７コマンドでコンパイルしてからバッチ処理すれば良い。ただし、ｆ７７コマンドのオプションは’-P　multi’，バッチジョブスクリプト中には’setenv　Ｆ＿RSVTASK　４’を指定する。
　場合によっては、これで並列化のための変更は完了という事もあるが、実際には、さらに変更が必要になると予想される。(第3段階)。それは、必要なデータをどのように割り付けるかという問題が残っているからである。特にデータを各CPUで共有するかどうかが重要な点である。例えば、固有値計算の場合、行列要素は（固有値計算が終われば必要でないならば）各CPUだけで有効であれば良いが、計算で得られた固有値は最終的にすべてのCPUで共有されるデータに通常格納されるであろう。このようなデータの区別はGLOBAL COMMON文やLOCAL COMMON文で宣言しておく必要がある。また、共有のデータへの格納のための文も当然必要となるので、さらなる変更が通常必要となる。この第３段階の変更が実際にはかなり重要で慎重に行なう必要があると考えられる。さて、共有データは通常GLOBAL COMMON文によって宣言するが、サブルーチンPTFORKの引数は宣言がなくても共有データとみなされる。第１段階でサブルーチンEIGENの引数が番号を指定する i だけになるように指定したのは、この点を配慮したからであり、共有データであれば引数として追加しても差し支えない（すなわち、以下のように変更しても良い）。

４　ｃａｌｌ　ｐｔｆｏｒｋ（ｉｔｖａｒ（ｉ），ｉｔｐｒａｍ（ｉ），ｅｉｇｅｎ，ｎｏ（ｉ），ａａ，ｂｂ，．．．
ｃａｌｌ　ｅｉｇｅｎ（ｎｏ（１），ａａ，ｂｂ，．．．）

　それでは、並列化の効果をＡＮＡＬＹＺＥＲで調べてみよう。ＡＮＡＬＹＺＥＲを起動するには、以下のコマンドをバッチジョブスクリプト中に含めて、バッチ処理すれば良い（ただし、最適化プリプロセッサを事前に起動させて、その結果をprog．fに保存しておく必要がある）。

ｆａｎｐ　－ｔｍ　－ｍｕｌｔｉ　ｐｒｏ．ｆ　/ｕｓｒ/ｌｉｂ/ａｓｌ　　＜　ｉｎｐｕｔ．ｄａｔ　＞　ｏｕｔｐｕｔ．ｄａｔ

　まず、マイクロタスク機能（自動並列化）を使った並列化から見てみよう。同じ計算（行列に時限は1600）をASLライブラリ中のDCSMAAおよびQCSMAA(並列版)サブルーチンを用いて行うと、計算の課金対称CPU時間（単位は秒）は、それぞれ、20.1および9.2となる。表１はANALYZER出力中の重要な数値を示したもので、上半分がDCSMAA，下半分がQCSMAAを用いた計算の結果である。それぞれの場合について、プログラム単位毎のCPU時間（正確には排他的CPU時間）、MFLOPS値、ベクトル演算率、平均ベクトル長、マイクロタスクレジデンス時間の全体に対する割合が示されている。全体のマイクロタスクレジデント時間は性能のMFLOPS値から見積もられるので、割合の後に時間の推定値も示してある。この割合が表示されていないプログラム単位ではCPU時間がマイクロタスクレジデント時間となる。ただし、CPU時間が全体の0.5％未満のプログラム単位については省いている。また、自動並列化を使っているので、プログラム単位が分割され、それらは元プログラム単位名＄１などで表されている。

表１：ANALYZERの結果１（マイクロタスク機能による並列化）
―――――――――――――――――――――――――――――――――――
PROG．UNIT	EXCLUSIVE CPU TIME	(　 　％)	MFLOPS	V．OP． RATIO	AVER． V．LEN	MICRO RES%
―――――――――――――――――――――――――――――――――――
HMLEVV TRANSEV$1 STRCHK STRCHK$1 EQIPAT$	17.321 6.738 0.367 0.312 0.273	（ 　68.9） （ 　26.8） （　　1.5） （　　1.2） （　　1.1）	1798.8 1945.3 0.9 427.1 47.8	99.62 99.23 0.09 99.69 96.53	218.9 228.5 68.3 165.8 184.3	8.6 0.4 0.4	1.68 0.08 0.08
―――――――――――――――――――――――――――――――――――
TRANSEV$1 HMLEVV STRCHK STRCHK$1 EQIPAT$	6.738 6.526 0.353 0.311 0.272	（　 47.2） （　 45.7） （　　2.5） （　　2.2） （　　1.9）	1945.4 1195.8 1.0 427.4 47.9	99.62 99.14 0.08 99.69 96.56	228.5 172.8 67.9 165.8 184.3	19.3 0.9 0.8	1.68 0.08 0.07
―――――――――――――――――――――――――――――――――――

表で、プログラム単位HMLEVVがDCSMAAあるいはQCSMAAサブルーチンに対応していると考えて良い。DCSMAAおよびQCSMAAサブルーチンでのCPU時間は、それぞれ、17.3および6.5程度、そのMFLOPS値は、両者のCPU時間の差(17.3－6.5＝10.8)がほぼ課金対象CPU時間の差(20.1－9.2＝10.9)となっている。その他のプログラム単位については、表の上下で、ほとんど差はない。これらのことから、QCSMAAサブルーチンが十分に並列化の機能を果たしていることが分かる。ただし、MFLOPS値はDCSMAAサブルーチンの２/３程度に低下しており、実質的な並列化の効果としては、４*２/３＝８/３程度である（17.3/6.5－8/3）。
　次に、マクロタスク機能を使った並列化を見てみよう。この場合は、DCSMAAサブルーチンを用いており、その課金対象CPU時間は10.1となる。この値は同じくDCSMAAを用いたマイクロタスク機能による並列化の計算での値20.1のほぼ半分である。出力例としてANALYZER出力およびプログラム情報の一部を最後に示しておくが、そのANALYZER出力はマイクロタスク機能による並列化の場合とは少し異なっている。特に注意しなければならないのはCPU時間であり、課金対称CPU時間との対応が異なっている。表２にANALYZER出力中の重要な数値を示すが、表１と比べる場合には、CPU時間はおおむねその１/４の値をとるべきであろうと思われる。

表２：ANALYZERの結果２（マクロタスク機能による並列化）
―――――――――――――――――――――――――――――――――
PROG．UNIT	EXCLUSIVE CPU TIME	(　 　％)	MFLOPS	V．OP． RATIO	AVER． V．LEN
―――――――――――――――――――――――――――――――――
HMLEVV TRANSEV EIGEN EQIPAT STRCHK	17.338 8.866 6.114 0.311 0.273	（ 　52.7） （ 　27.0） （　 18.6） （　　0.9） （　　0.7）	1797.1 1478.3 20.8 49.0 0.7	99.62 99.36 73.54 96.73 0.12	218.9 228.5 14.7 185.9 82.2
―――――――――――――――――――――――――――――――――

表より、プログラム単位HMLEVVのMFLOPS値が表１のDCSMAAサブルーチンを用いた場合の値にほぼ等しくなっているのが分かる。しかし、プログラム単位TRANSEVのMFLOPS値は表１の値より低下している。また、並列化だけのためのプログラム単位EIGENにかなりのCPU時間が費やされている。マクロタスク機能による並列化には、高いMFLOPS値を持つDCSMAAサブルーチンを用いるという利点もあるが、その反面、並列化に必要なCPU時間の増加などの欠点も生じてくる。
　さて、ここまでの計算だけを見ると、課金対象CPU時間はQCSMAAを用いた計算の9.2が最小である。しかし、出力例で気になる点がある。それは、プログラム情報で、４台以上で実行した時間が約4.0と、１台以上で実行した時間などと比べると、極端に小さいことである。これは、固有値計算に要する時間が並列によって異なることが原因である。この計算では４個の行列の固有値を計算しているが、実は、そのうち１個の行列の固有値計算が他に比べてかなり容易であったのである。そこで、同じような行列４個の固有値計算を行うと、QCSMAAを用いた計算では11.3、マクロタスク機能による並列化による計算のほうが優るようになった。一般的には、この結果のほうが正しいといえる。
　結論として、マクロタスク機能による並列化はかなり有効である。並列版のサブルーチンが用意されている場合でも、マクロタスク機能による並列化のほうが優っている。またCPU数を増加させても、その効率が極端に低下することはあまりないと予想される。