Tutorials:MonteCarloHOWTO/ja

概要

ALPSアプリケーションでは、すでにいくつかのモンテカルロ（MC）プログラムが入手可能です。ALPSの主な目的は、開発者が自身のMCプログラムを簡単に、短時間で開発可能なように手助けをすることです。ALPSアプリケーションには、すべてのMCアルゴリズムで実行されるいくつかの普遍的な機能があります。もう再開発の必要はありません。この章では、ALPSライブラリが提案するシンプルなフレームワークを説明し、古典モンテカルロ、量子モンテカルロ法についていくつかの例で演習をおこないます。

ALPSソースのexample/schedulerディレクトリにイジングモデルのサンプルがあります。このサンプル計算例を用いて、以下で解説します。

　ALPSライブラリを使用する利点　

• 　エラーバーと自己相関時間の自動計算
• 　自動並列化
• 　抽出ツールを使用した計算結果の自動出力（XML形式）
• 　符号の自動制御（量子モンテカルロシミュレーションの符号問題）
• 　簡単なチェックポイント、再計算機能
• 　入力パラメータの容易な制御
• 　新しい測定を容易に追加
• 　乱数の容易な取扱い
• 　格子の容易な取扱い
• 　量子ハミルトニアンの容易な取扱い（量子モンテカルロ計算）
• ...

わずかなプログラミング/移植によって、これらの利点が得られます。

ALPS MCプログラム使用の第一歩

このチュートリアルは未完成です。コメント、ご提案があれば"comp-phys-alps-users@lists.comp-phys.org"までお願いします。

First of all, using the ALPS libraries implies that you use C++ as the programming language. From now on, let us assume that you want to program a new MC application and that you already know what are going to be your internal data structure and the algorithm you will use for a single Monte Carlo step. まず最初に、ALPSライブラリを使用するためにC＋＋言語をプログラム言語としてください。今から、新しいMCアプリケーションの開発をおこない、シングルモンテカルロステップで使用される内部データ構造やアルゴリズムは分かっているものと仮定します。

ALPSライブラリを使うために、新しいC＋＋クラスが必要です。このクラスは内部ALPSクラスとは分離されています。このクラスをMyMonteCarloと呼び、MyMonteCarlo.hppヘッダーファイルを生成してください。

 #ifndef MYMC_HPP
 #define MYMC_HPP
 #include <alps/scheduler.h>
 #include <alps/alea.h>
 #include <alps/scheduler/montecarlo.h>
 #include <alps/osiris.h>
 #include <alps/osiris/dump.h>
 #include <alps/expression.h>
 using namespace std;
 class MyMonteCarlo : public alps::scheduler::MCRun{
 public :
   MyMonteCarlo(const alps::ProcessList &,const alps::Parameters &,int);
   static void print_copyright(std::ostream &);
   void save(alps::ODump &) const;
   void load(alps::IDump &);
   void dostep();
   bool is_thermalized() const;
   double work_done() const;
 private :
 // your own internal data here ...
 };
 #endif

各ラインが何を意味するか説明していきます。

• 　最初に、ALPSのヘッダーファイルをインクルードする必要があります。（他のALPS関数を利用する場合は、さらに追加する必要あり）
• 　ALPSクラスから派生させることによってMyMonteCarloクラスを定義します。上記のヘッダーファイルは、最も簡単なクラスであるMCRun ALPSクラスを使用することを意味しています。もしALPSの格子機能を使用したい場合は、次の用に置き換えてください。

 class MyMonteCarlo : public alps::scheduler::LatticeMCRun<>{

ヘッダーの追加もおこなってください。

 #include <alps/lattice.h>

もし、さらにモデルライブラリ（格子量子モデル）を使用する場合は、次の用に置き換えてください。

 class MyMonteCarlo : public alps::scheduler::LatticeModelMCRun<>{

• 　publicにある、MyMonteCarlo(const alps::ProcessList&,const alps::Parameters&,int)はコンストラクターです。そして、すべてのパラメータの初期化が必要でしょう。
• 　print_copyright(std::ostream&)は、各シミュレーションの最初の部分をアウトプットするかどうか設定する関数です。
• 　save,loadは、各チェックポイントの後で再スタートするために、何がディスク上にsave(load)する必要があるか決定する関数です。
• 　dostep関数はMC関数のメインです。あらゆるMDステップで実行されます。
• 　is_thermalizedは熱化のパートのシミュレーションがいつ終わったかALPSライブラリに教える関数です。
• 　work_doneは、シミュレーションの何パーセントが実行済みかALPSライブラリに教える関数です。　

最後に、これらの関数とプログラムを正しく結びつけるプログラム、MyMonteCarlo.cppを準備します。

ALPS MCプログラムのビルド

copyright関数

まず、MyMonteCarlo.cppの中で最も簡単なprint_copyright()関数から始めます。

 #include "MyMonteCarlo.hpp";
 /************************************ ALPS functions **********************************************/
 // Copyright statement
 void MyMonteCarlo::print_copyright(std::ostream & out)
 {
   out << "My own ALPS Monte Carlo program v. 0.1\n"
       << "  copyright (c) 2006 by Myself,\n"
       << "  available from the author on request\n\n";
 }

モンテカルロステップの制御

次に、MCステップについて説明します。典型的な例は次の通りです。まず、熱化のためステップ数を決め、続いて、測定領域のステップ数を決めます。普通、各計測の間にモンテカルロステップを実行します。よって、内部データ構造で、次の変数(in MyMonteCarlo.hpp)を設定することで実行可能です。

 private :
  // your own internal data here ...
  int Nb_Steps; 
  int Nb_Thermalisation_Steps; 
  int Each_Measurement;
  int Steps_Done_Total; 
  int Measurements_Done;
  void do_update();
  void do_measurements();
  // the rest of your own internal data here ...
 };

Nb_Thermalisation_Stepsは熱化ステップ数、Nb_Stepsは熱化後のMCステップ数です。Each_Measuremenは各測定のステップ数をあらわします。これらの３つの数字はコンストラクタで初期化されます。Steps_Done_TotalはMCが終わったステップの番号を保存します。最後に、Measurements_Done willは各測定の計算し終わった途中のステップ数を保存します。The do_update()とdo_measurements()関数はそれぞれ各一回のMCステップで実行されます。

これらすべての定義はMyMonteCarlo.cppでクラスメンバー関数のシンプルな定義となります。

 bool MyMonteCarlo::is_thermalized() const
   {  return (Steps_Done_Total >= Nb_Thermalisation_Steps); }

この間数は、現在のMCステップ数が熱化ステップのトータル番号より大きいとき、１を返し、それ以外は０を返します。

 double MyMonteCarlo::work_done() const
 { return (is_thermalized() ? (Steps_Done_Total-Nb_Thermalisation_Steps)/double(Nb_Steps) :0.); }

シミュレーションが熱化しなければ０を返し、測定ステップがどれくらい実行されているかパーセンテージを与えます。最後に、dostep()関数は次のようになります。

 void MyMonteCarlo::dostep()
 { do_update(); // you'll have to define what this function does later
   ++Steps_Done_Total; // increment the number of steps done
   if (is_thermalized()) // do measurements only if simulation thermalized
     { if (++Measurements_Done==Each_Measurement) // do a measurement every Each_Measurement
       { Measurements_Done=0;
       do_measurements(); // you'll have to define what this function does later
       }
    }
 }

いくつかのポイントで、do_update()によっておこなわれる、MCステップ中に何をおこなうのか？ちゃんと調べなければならない。計測は例にあるdo_measurements()関数を参照してください。

セーブ、ロード関数

ALPSには計算中の内部データをディスクにチェックポイントを生成する機能があります。チェックポイントほ利用して計算を再スタートするために、整数データと倍精度の配列が必要になります。

 private :
  // the rest of your own internal data here ...
  int Number_of_Spins; 
  int MyOwnVariable;
  std::vector<double> SpinArray;
  ...
 };

設定は簡単です。

 void MyMonteCarlo::save(alps::ODump& dump) const
   { dump <<  Number_of_Spins << MyOwnVariable << SpinArray;}

load関数は簡単に推測できます。

 void MyMonteCarlo::load(alps::IDump& dump) const
   { dump >>  Number_of_Spins >> MyOwnVariable >> SpinArray;}

ここで、ダンプの方法を説明します。一般的な型（int,double,boolなど）や一般的なコンテナ（vectorやset）が対応可能です。

次のような内部構造を持っているとすると、

 struct Vertex
   { int vertex_type;
    std::vector<double> coordinates;
    int SomeOtherVariable; }

モンテカルロクラスでインスタンス化し保存したいと考えると、

 private :
    // the rest of your own internal data here ...
    Vertex MyVertex; 
    ...
   };

構造の中で、何を保存、ロードするのかをALPSに指示しなくてはいけません。

 alps::ODump& operator<<(alps::ODump& dump, const Vertex& v)
 { return dump << v.vertex_type << v.coordinates; }

 alps::IDump& operator>>(alps::IDump& dump, Vertex& v)
 { return dump >> v.vertex_type >> v.coordinates;}

そして、ALPSのsave,load関数のMyVertexを追加することが容易にできます。

 void MyMonteCarlo::save(alps::ODump& dump) const
   { dump <<  Number_of_Spins << MyOwnVariable << SpinArray << MyVertex;}

もし、上記のようにモンテカルロステップ数を設定したならば、save/load関数は次の用に加えることができます。

 void MyMonteCarlo::save(alps::ODump& dump) const
   { dump <<  Number_of_Spins << MyOwnVariable << SpinArray << MyVertex;
     dump <<  Nb_Steps << Measurements_Done; }

do_measurements()関数

物理量の計測の実行方法です。次の例を参照してください。

 void MyMonteCarlo::do_measurements()
   { double Energy; std::valarray<double> Correl(L);
     // do the measurements in your code (update the Energy and Correl variable) ...
     ...

     // give them to ALPS
     measurements["Energy"] << Energy;
     measurements["Spin Correlations"] << Correl;
   }

わかりますか？驚くかもしれませんが、ALPSは "Energy"や"Spin Correlations"をどうやって知るのか？また、スカラー値（エネルギーなど）とベクター値（相関）をどうやって区別するのか？この２つの疑問は、クラスのコンストラクタで以下の用に対処されます。

他では、なぜベクター計測値にstd::vectorではなく、std::vectorを使用しているのか？これはALPS内部の問題に起因しています。ALPSは各測定値で常に同じサイズである必要があるからです。