scan2d

調和拘束と MLIP 緩和により、2 距離 (d₁, d₂) のグリッドスキャンを行い、(d₁, d₂) 上の 2D ポテンシャル面を得ます。TS 領域の特定や、MEP 精密化前の反応ランドスケープの可視化に使います。scan2d--max-step-size に基づいて両軸の線形グリッドを作成し、各格子点を対応する拘束付きで緩和して、可視化用にバイアスなしの MLIP エネルギーを記録します。入力は 1 つの構造 + -s/--scan-lists scan2d.yaml(推奨)、またはちょうど 2 つの 4 要素タプルを含む -s/--scan-lists単一 インラインリテラルです。デフォルトのバックエンドは UMA で、-b/--backend で他のバックエンドも選択できます。L-BFGS の代わりに RFOptimizer を使う場合は --opt-mode hess を指定してください。

XYZ/GJF 入力では、--ref-pdb で参照 PDB トポロジーを指定すると、XYZ 座標を保持したまま PDB/GJF へのフォーマット対応変換が可能になります。

実行例

YAML スペックファイルを使った最小実行:

pdb2reaction scan2d -i input.pdb -q 0 -s scan2d.yaml -o ./result_scan2d/

コマンド形式:

pdb2reaction scan2d -i INPUT.{pdb|xyz|trj|...} [-q CHARGE] [-l, --ligand-charge <number|'RES:Q,...'>] [-m MULT] \
 [-b/--backend uma|orb|mace|aimnet2] [--solvent SOLVENT] [--solvent-model alpb|cpcmx] \
 [-s/--scan-lists scan2d.yaml | '[(i,j,lowÅ,highÅ), (i,j,lowÅ,highÅ)]'] [options] \
 [--convert-files/--no-convert-files] [--ref-pdb FILE]

推奨: YAML/JSON スペックファイル:

cat > scan2d.yaml << 'YAML'
one_based: true
pairs:
 - ["TYR,285,CA", "SAM,309,C10", 1.30, 3.10]
 - ["TYR,285,CB", "SAM,309,C11", 1.20, 3.20]
YAML
pdb2reaction scan2d -i input.pdb -q 0 -s scan2d.yaml

代替: インライン Python リテラル:

pdb2reaction scan2d -i input.pdb -q 0 \
 -s '[("TYR,285,CA","SAM,309,C10",1.30,3.10),("TYR,285,CB","SAM,309,C11",1.20,3.20)]'

LBFGS、内側軌跡ダンプ、Plotly 出力:

pdb2reaction scan2d -i input.pdb -q 0 \
 -s '[("TYR,285,CA","SAM,309,C10",1.30,3.10),("TYR,285,CB","SAM,309,C11",1.20,3.20)]' \
 --max-step-size 0.20 --dump -o ./result_scan2d/ --opt-mode grad \
 --preopt --baseline min

スキャンリスト仕様

scan2d はちょうど 2 つの4 要素タプル (i, j, low_Å, high_Å) を受け付けます(YAML/JSON では pairs キー、インラインでは単一リテラル)。scan と異なり、リテラルは 1 つだけを受け付けます(複数ステージは非対応)。

YAML/JSON ファイル書式、インライン Python リテラル構文、原子セレクタ、クォート規則については CLI 規約: スキャンリスト仕様 を参照してください。

処理の流れ

  1. geom_loader で入力構造をロードし、電荷とスピンを解決します。--preopt の場合は無バイアスの事前最適化を実行します。-q が省略され --ligand-charge がある場合、構造は酵素–基質複合体として扱われ、PDB 入力(または --ref-pdb 付き XYZ/GJF)では extract.py の電荷サマリーから総電荷を導出します。事前最適化構造は grid/preopt_iDDD_jDDD.*DDD = round(d × 100) Å)に保存され、surface.csv には i = j = -1 のエントリとしてバイアスなしエネルギーが記録されます。

  2. -s/--scan-lists(YAML/JSON ファイルパスまたはインライン Python リテラル)を 2 つの4 要素タプルに解析し、インデックスを正規化します(デフォルトは 1 始まり)。PDB 入力では、各エントリに整数インデックスまたは 'TYR,285,CA' のようなセレクタ文字列を指定できます。区切りは空白・カンマ・スラッシュ・バッククォート・バックスラッシュのいずれも可で、トークン順序は任意です(フォールバックは resname, resseq, atom を想定)。線形グリッドは ceil(|high low| / h) + 1 点(両端を含む)で構成します(h = --max-step-size)。長さ 0 の範囲は 1 点に縮退します。その後、各軸は事前最適化構造に最も近い距離が i = 0 / j = 0 になるよう並べ替えられます。

  3. 外側ループで d1[i](近い順)を走査します。各値で d₁ 拘束のみを適用して緩和し、その構造をスナップショットとして保存します。次に内側ループで d2[j] を走査し、d₁ と d₂ の両拘束を適用して、最も近い既収束構造から緩和を開始します。

  4. (i, j) について、<out-dir>/grid/point_iDDD_jDDD.xyzDDD = round(d × 100) Å。例 d1=1.30 Å, d2=3.10 Åpoint_i130_j310.xyz)に構造を保存し、バイアス収束の可否を記録し、バイアスを除去した MLIP エネルギーを評価します。--dump の場合、外側ループごとの内側軌跡が inner_path_d1_###_trj.xyz### は外側ステップ index)として保存されます。

  5. すべての点を走査したら、i,j,d1_A,d2_A,energy_hartree,bias_converged,energy_kcal,d1_label,d2_label の列を持つ <out-dir>/surface.csv を作成し、--baseline {min|first} で kcal/mol の基準をシフトします。--baseline first では基準点が (low₁, low₂) ではなく再並べ替え後の最初の格子点(i = j = 0)になります。scan2d_map.png(2D コンター)と scan2d_landscape.html(3D サーフェス)を <out-dir>/ に生成します。--zmin/--zmax でカラースケールを固定できます。

出力

実行後は surface.csvgrid/ 配下の各点の構造、scan2d_map.png / scan2d_landscape.html のプロットを確認してください。

  • result_scan2d/surface.csv

  • result_scan2d/grid/point_iDDD_jDDD.xyz (DDD = round(d × 100) Å。例: d1=1.30 Å, d2=3.10 Åpoint_i130_j310.xyz)

  • result_scan2d/scan2d_map.pngresult_scan2d/scan2d_landscape.html

out_dir/ (デフォルト:./result_scan2d/)
├─ surface.csv # 構造化グリッド表
├─ scan2d_map.png # 2D コンター(Kaleido 必須; PNG 出力に失敗すると実行が停止)
├─ scan2d_landscape.html # 3D サーフェス可視化
├─ grid/point_iDDD_jDDD.xyz # DDD = round(d × 100) Å (例: d1=1.30 Å, d2=3.10 Å → point_i130_j310.xyz)
├─ grid/point_iDDD_jDDD.pdb # 変換有効時に対応する PDB
├─ grid/point_iDDD_jDDD.gjf # テンプレートがある場合に対応する Gaussian
├─ grid/preopt_iDDD_jDDD.xyz # 事前最適化構造(--preopt が True の場合)、DDD = round(d × 100)
├─ grid/preopt_iDDD_jDDD.pdb # 変換有効時に対応する PDB
├─ grid/preopt_iDDD_jDDD.gjf # テンプレートがある場合に対応する Gaussian
└─ grid/inner_path_d1_###_trj.xyz # --dump の場合のみ(### は外側ステップ index、PDB 入力時は.pdb にも変換)

CLI オプション

オプション

説明

デフォルト

-i, --input PATH

geom_loader が受け入れる構造ファイル

必須

-q, --charge INT

総電荷(CLI > テンプレート/--ligand-charge)。両方指定時は -q が優先

テンプレート/導出がない場合は必須

-l, --ligand-charge TEXT

単一の整数(例: -1)でリガンド総電荷を指定するか、残基別マッピング(例: GPP:-3,SAM:1)で PDB 残基電荷から全系の電荷を導出。-q 省略時に使用(PDB 入力、または --ref-pdb 付き XYZ/GJF)

None

--workers, --workers-per-node

MLIP 予測器の並列度(workers > 1 で解析Hessian無効; UMA バックエンドのみ; workers_per_node は並列予測器へ転送)。診断上の注意は workers > 1 は解析Hessianを無効化する(UMA バックエンド) を参照

1, 1

-m, --multiplicity INT

スピン多重度 2S+1。.gjf テンプレートがあれば継承し、未指定時は 1

.gjf テンプレート値または 1

-s, --scan-lists TEXT

スキャンターゲット: YAML/JSON スペックファイルパス(推奨)または 単一のインライン Python リテラルで 2 つの4 要素タプル (i,j,lowÅ,highÅ) を指定。i/j は整数インデックスまたは PDB セレクタ('TYR,285,CA'

必須

--one-based/--zero-based

(i, j) のインデックスを 1 始まり/0 始まりとして解釈

True

--print-parsed/--no-print-parsed

-s/--scan-lists 解釈後のペア情報を表示

False

--max-step-size FLOAT

各距離の 1 増分あたりの最大変化量(Å)。グリッド密度を決定

0.20

--bias-k FLOAT

調和バイアス強度 k(eV·Å⁻²)

300

--relax-max-cycles INT

各バイアス緩和の最大最適化サイクル数。YAML で opt.max_cycles が指定されていない場合に使用

10000

--opt-mode TEXT

grad → L-BFGS、hess → RFOptimizer

grad

--freeze-links/--no-freeze-links

PDB 入力時にキャップ水素の親原子を凍結

True

--freeze-atoms TEXT

凍結する原子の 1 始まりインデックスをカンマ区切りで明示的に指定(例: '1,3,5')。--freeze-links と併用可、任意の入力形式に適用

None

--dump/--no-dump

外側ループごとの inner_path_d1_###_trj.xyz を保存

False

--convert-files/--no-convert-files

PDB/Gaussian 入力で XYZ/TRJ → PDB/GJF 変換を切り替え

True

--ref-pdb FILE

XYZ/GJF 入力時の参照 PDB トポロジー(XYZ 座標を保持)

None

-o, --out-dir TEXT

出力ディレクトリ

./result_scan2d/

--thresh TEXT

収束プリセットの上書き(gau_loose, gau, gau_tight, gau_vtight, baker, never

baker

--config FILE

ベース YAML 設定ファイル(最初に適用)

None

-b, --backend {uma,orb,mace,aimnet2}

MLIP バックエンド

uma

--solvent TEXT

xTB 暗黙溶媒(例: water)。none で無効化

none

--solvent-model {alpb,cpcmx}

xTB 溶媒モデル

alpb

--preopt/--no-preopt

スキャン前に無バイアス最適化を実行

False

--baseline {min,first}

kcal/mol の基準をグローバル最小値または最初の格子点に設定

min

--zmin FLOAT, --zmax FLOAT

カラースケールの下限/上限(kcal/mol)

自動

--out-json/--no-out-json

out_dir に機械可読な result.json を書き出す。スキーマは JSON 出力スキーマ を参照

False

YAML 設定

geom:
 coord_type: cart # coordinate type: cartesian vs dlc internals
 freeze_atoms: [] # 1-based frozen atoms merged with CLI/link detection
calc:
 charge: 0 # total charge (CLI/template override)
 spin: 1 # spin multiplicity 2S+1
 model: uma-s-1p1 # uma-s-1p1 | uma-m-1p1
 device: auto # MLIP device selection
opt:
 thresh: baker # convergence preset (default: baker)
 max_cycles: 10000 # optimizer cycle cap
 dump: false # optimizer dumps (scan trajectories are controlled by --dump)
 out_dir: ./result_scan2d/ # output directory
lbfgs:
 max_step: 0.3 # maximum step length
 out_dir: ./result_scan2d/ # LBFGS-specific output directory
rfo:
 trust_radius: 0.10 # trust-region radius
 out_dir: ./result_scan2d/ # RFO-specific output directory
bias:
 k: 300.0 # harmonic bias strength (eV·Å⁻²)

共有 YAML セクション

  • geom, calc, opt, lbfgs, rfo: YAML リファレンス と同じキーを使用します。opt.dump は YAML で設定可能ですが、スキャン軌跡の出力は --dump で制御します。

opt の詳細は YAML リファレンス を参照してください。

注意事項

  • 計算エンジンは MLIP バックエンド(デフォルト: UMA、-b/--backend で切替可能)で、1D スキャンと同じ HarmonicBiasCalculator を再利用します。

  • Å 単位の制限値は内部で Bohr に変換され、L-BFGS ステップや RFO 信頼半径の制御に使われます。最適化の一時ファイルはテンポラリディレクトリに配置されます。

  • バイアスはエネルギー記録前に除去されるため、surface.csv を下流のフィッティングや可視化スクリプトにそのまま利用できます。

  • --freeze-links はユーザー指定の freeze_atoms にキャップ水素親原子をマージし、抽出された活性部位モデルの境界を固定します。

  • --relax-max-cycles明示的に指定され、かつ YAML で opt.max_cycles が設定されていない場合にのみ適用されます(デフォルト 10000)。

関連項目