チュートリアル08:water finger の解析と2次元自由エネルギー面の計算

本章ではwater finger の解析と2次元自由エネルギー面の計算について説明する。

計算する系の作成

　チュートリアル07で用いた系と同様の系を考える。Clイオンと水層の重心とのz方向の距離を一般化座標zとし、water finger を射影した座標を一般化座標wとする。 water finger を解析するための初期構造は、イオンが界面油相側に存在する方が扱いやすい。そこでまず、Cl^-を水重心から20Å離れた位置に置く。系はチュートリアル07で用いたものを使えばよい。

ネームリストの作成①（系の平衡化）

water finger座標wの設定

　チュートリアル07で用いたネームリストにwater finger座標wを加える。

=== input.nml ================================================================================
&FreeFlex
  (略)
/

&spcoord
  condition = "WATFING 878 877 -0.4d-10"
/
===============================================================================================

バイアスの設定

　二次元のバイアスを用意する。バイアスの置き方は以下に示す。バイアス番号さえ対応していれば，図のようなピラミッド型ではなく縦に並べてもよい。

=== input.nml =================================================================================

  (略)
 
&biasdata
PARABOLIC 1 1 2 0.5d-10 0.1000000000E-08  1
PARABOLIC 1 1 2 0.5d-10 0.1100000000E-08  2
PARABOLIC 1 1 2 0.5d-10 0.1200000000E-08  3
PARABOLIC 1 1 2 0.5d-10 0.1300000000E-08  4
PARABOLIC 1 1 2 0.5d-10 0.1400000000E-08  5
PARABOLIC 1 1 2 0.5d-10 0.1500000000E-08  6
PARABOLIC 1 1 2 0.5d-10 0.1600000000E-08  7
PARABOLIC 1 1 2 0.5d-10 0.1700000000E-08  8
PARABOLIC 1 1 2 0.5d-10 0.1800000000E-08  9  25
PARABOLIC 1 1 2 0.5d-10 0.1900000000E-08  10 26 41
PARABOLIC 1 1 2 0.5d-10 0.2000000000E-08  11 27 42 56
PARABOLIC 1 1 2 0.5d-10 0.2100000000E-08  12 28 43 57 70
PARABOLIC 1 1 2 0.5d-10 0.2200000000E-08  13 29 44 58 71 83
PARABOLIC 1 1 2 0.5d-10 0.2300000000E-08  14 30 45 59 72 84 95
PARABOLIC 1 1 2 0.5d-10 0.2400000000E-08  15 31 46 60 73 85 96  106
PARABOLIC 1 1 2 0.5d-10 0.2500000000E-08  16 32 47 61 74 86 97  107 116
PARABOLIC 1 1 2 0.5d-10 0.2600000000E-08  17 33 48 62 75 87 98  108 117 125
PARABOLIC 1 1 2 0.5d-10 0.2700000000E-08  18 34 49 63 76 88 99  109 118 126 133
PARABOLIC 1 1 2 0.5d-10 0.2800000000E-08  19 35 50 64 77 89 100 110 119 127 134 140
PARABOLIC 1 1 2 0.5d-10 0.2900000000E-08  20 36 51 65 78 90 101 111 120 128 135 141 146
PARABOLIC 1 1 2 0.5d-10 0.3000000000E-08  21 37 52 66 79 91 102 112 121 129 136 142 147 151
PARABOLIC 1 1 2 0.5d-10 0.3100000000E-08  22 38 53 67 80 92 103 113 122 130 137 143 148 152 155
PARABOLIC 1 1 2 0.5d-10 0.3200000000E-08  23 39 54 68 81 93 104 114 123 131 138 144 149 153 156 158
PARABOLIC 1 1 2 0.5d-10 0.3300000000E-08  24 40 55 69 82 94 105 115 124 132 139 145 150 154 157 159 160
PARABOLIC 1 2 2 0.5d-10 0.3000000000E-09  25  26  27  28  29  30  31  32  33  34  35  36 37 38 39 40
PARABOLIC 1 2 2 0.5d-10 0.4000000000E-09  41  42  43  44  45  46  47  48  49  50  51  52 53 54 55
PARABOLIC 1 2 2 0.5d-10 0.5000000000E-09  56  57  58  59  60  61  62  63  64  65  66  67 68 69
PARABOLIC 1 2 2 0.5d-10 0.6000000000E-09  70  71  72  73  74  75  76  77  78  79  80  81 82
PARABOLIC 1 2 2 0.5d-10 0.7000000000E-09  83  84  85  86  87  88  89  90  91  92  93  94
PARABOLIC 1 2 2 0.5d-10 0.8000000000E-09  95  96  97  98  99  100 101 102 103 104 105
PARABOLIC 1 2 2 0.5d-10 0.9000000000E-09  106 107 108 109 110 111 112 113 114 115
PARABOLIC 1 2 2 0.5d-10 0.1000000000E-08  116 117 118 119 120 121 122 123 124
PARABOLIC 1 2 2 0.5d-10 0.1100000000E-08  125 126 127 128 129 130 131 132
PARABOLIC 1 2 2 0.5d-10 0.1200000000E-08  133 134 135 136 137 138 139
PARABOLIC 1 2 2 0.5d-10 0.1300000000E-08  140 141 142 143 144 145
PARABOLIC 1 2 2 0.5d-10 0.1400000000E-08  146 147 148 149 150
PARABOLIC 1 2 2 0.5d-10 0.1500000000E-08  151 152 153 154
PARABOLIC 1 2 2 0.5d-10 0.1600000000E-08  155 156 157
PARABOLIC 1 2 2 0.5d-10 0.1700000000E-08  158 159
PARABOLIC 1 2 2 0.5d-10 0.1800000000E-08  160
DUMMY     1 2 1 0.5d-10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
&end
===============================================================================================

計算の実行①（系の平衡化）

　チュートリアル07と同様にして並列計算を行う。100ps ほど平衡化を行えばよい。

> mpiexec –machinefile machinefile –n 160 ./FreeFlex.exe input.nml &

ネームリストの作成② （REUS計算）

REUSの設定

　二次元のREUSを行うためにネームリストを編集する。

=== input2.nml =================================================================================

  (略)
 
&REUS
  calculate_reus = .TRUE.
  dim = 2
  max_iter = 100
  exchange_interval = 100
  bias_set_id = 1
  dim_q = "1 1"
  pos_PARA = "2 2"
  spcoord_id = "1 2"
/
===============================================================================================

データ出力の設定

　後でWHAM (重み付きヒストグラム解析法)を行うためのデータを出力する。

=== input2.nml =================================================================================

  (略)
 
&MONITOR
  interval = 1
  file = "z.output"
  condition = "REUS_qU 1 2"
  separator = "$"
/
===============================================================================================

計算の実行②（REUS計算）

　チュートリアル02と同様にして、100psほど並列計算を行う。

> mpiexec –machinefile machinefile –n 160 ./FreeFlex.exe input2.nml &

WHAMの実行

　アンブレラサンプリングを行ったデータを解析するにはWHAMを用いる必要がある。WHAMを実行するためには、先のネームリスト"input.nml"の&MONITORを&WHAMSETTINGと&WHAM_INPUTに変えたネームリストを作らなければならない。

=== wham.nml ==================================================================================

  (略)
 
&WHAMSETTING
  wham_mode = .TRUE.
  mpi_size = 160
/

&WHAM_INPUT
  wham_mode = .TRUE.
  reus_mode = .TRUE.
  show = .TRUE.
  input_file_name = "z.output"
  free_energy_output_file_name = "free_energy.output"
  histogram_output_file_name = "histogram.out"
  nbin = 200
  ndir = 160
  ncore = 160
  nprof_column = 1
  data_column = "3 4"
  spcoord_ID = "1 2"
  bias_set_ID = 1
  interval = 1
  nskip = 14
  nread = 1000000
/
===============================================================================================

これを実行する。このときMPIは用いないことに注意する。

> ./FreeFlex.exe wham.nml &

計算の結果

　上記の計算を実行することによって以下のような一次元の自由エネルギー面と二次元の自由エネルギー面が得られる。以下の図において界面の位置はそれぞれ $z=0$ である。