技術論文
3Dプリンターによる分子キャビティーの構造解析
2016 年 15 巻 4 号 p. 97-104
詳細
2016 年 15 巻 4 号 p. 97-104
立体規則性高分子の単結晶において，分子鎖間に「分子キャビティー」と呼ばれる空孔が生成する系があり，その構造を制御することにより，特定の分子を分離・回収する高性能な分離膜としての利用が期待されている．分子キャビティーの構造を解明するために様々な研究が行われているが，複雑なキャビティーの構造を立体的・直観的に把握することは困難である．そこで，本研究では，3Dプリンターを用いてキャビティーの立体モデルを作成することで，構造の理解を深めることを目的とする．分子キャビティーを持つシンジオタクチックポリスチレン(s-PS)のδe型，ε型及びS-I型を対象として，分子動力学(MD)シミュレーションを行い，分子キャビティーの構造を再現した．得られたデータを3Dプリンターで利用可能な構造データに変換し，分子キャビティーの立体造形モデルを作成した．δe型について，b軸方向及びac軸方向につながるチャネルを再現したところ，気体分子が拡散しやすいac軸方向につながるチャネルの方が太くて短い様子が十分に再現されていた．一方ε型とS-I型については，それぞれ，円筒型とジグザグ状のチャネルが得られた．さらに，これらのモデルに気体分子を挿入した系について解析を行い，気体分子が跳躍運動をしている瞬間を3Dモデル化することにより，結晶構造と気体の拡散様式の関係を直観的にとらえることに成功した．
近年，地球表面の大気や海洋の平均温度が上昇する地球温暖化が問題となっており，二酸化炭素やメタンなどの温室効果ガスがその原因とされている．そこで，地球温暖化対策の一つとして期待されているのが，二酸化炭素回収・貯留技術(Carbon dioxide Capture and Storage; CCS)である．CCSとは，火力発電所や製鉄所のように二酸化炭素が多く排出されるような場所で，大量の燃焼ガスの中から二酸化炭素だけを回収し，地中に貯留または海中に注入する技術のことである．このように温室効果ガスのような低分子だけを効率的に分離・回収する技術として，本研究では膜分離法に注目する．
分離膜の分子設計において，選択性と透過量の両立が求められる．この2つの要件を満たすために高分子結晶を用いるアイデアが提案されている [1,2,3,4]．本来，高分子結晶は密度が高く，あまり間隙が空いていないものであるが，主鎖の立体規則性(タクティシティー)が制御された高分子の結晶において，分子鎖間に「一定の大きさと形状及び秩序的結合性を有する，ナノスケールの空孔」が形成される場合がある．このような空孔を「分子キャビティー」と呼び，結晶中の分子キャビティーの大きさ，形状，結合性を精密に制御することにより，特定の気体のみを吸着，あるいは，選択透過する，高機能な分離膜(スマートメンブレン) [1]として応用することができると考えられている．分子キャビティーの構造を解明するために，様々な研究が行われているが，複雑なキャビティーの構造を，立体的・直観的に把握することは困難である．
そこで本研究では，近年，複雑な造形を実現できるとして，工業・化学・医療などの分野で利用されている3Dプリンターを活用し，分子動力学(MD)シミュレーションと組み合わせることにより，分子キャビティーのモデルを作成し，キャビティー構造の理解を深める．また，気体分子がキャビティー間のチャネルをジャンプする瞬間のモデルを作成し，気体の拡散機構解明の一助とする．
著者(玉井)が開発した分子シミュレーションプログラムPAMPSを使用して，MDシミュレーション及び自由体積解析を行い，得られた構造データをもとに，3Dプリンターを用いて分子キャビティーの造形を行った．以下にその手順を示す．
2.1 分子動力学シミュレーション初期構造は，X線回折実験 [5,6]により決定された部分座標をもとに，空間群などの情報を用いて生成した．MDシミュレーションの基本セル中には，Table 1に示す数の結晶単位胞が含まれ，総原子数はどちらも9216 (高分子鎖数24本)である．三次元周期境界条件により，無限に広がった単結晶モデルとした．
| Form | Space group | Nunit | Natom |
| δe | P21/a[5] | 3×4×6 | 9216 |
| ε | Pbcn[6] | 3×2×6 | 9216 |
生成した初期構造に対して，最初に，原子間の過度な重なりを防ぐためにエネルギーの極小化を行った．50 psの平衡化シミュレーションの後，100 psのサンプリングシミュレーションを行った．平衡化シミュレーションでは，10 psごとにMaxwell分布によって各原子に速度をランダムに与えた．運動方程式は，時間刻み幅1 fsで数値積分した．相互作用計算のカットオフ半径は12.0 Åとし，クーロンポテンシャルに関しては，Ewald法を用いて無限遠まで考慮した．温度制御には能勢の方法，圧力制御にはパリネロ・ラーマン法を使用した．この方法では，セルの変形が許されており，安定な結晶構造を再現することができる．パラメータは汎用力場Amber (ff99)を合成高分子に拡張したもの [4,7]を用いた．設定温度は300 K，圧力は1.0 × 105 Paとした．
また，著者(玉井)の既報 [8]の手順に従い，ε型結晶に対してb軸方向の応力を印加することにより，構造転移を経て，サブナノスケールのキャビティー構造をもつS-I型結晶を得た．さらに，各結晶中にCO2またはN2を24分子挿入した系についても，同様の条件でシミュレーションを行い，それぞれ，1 ns のサンプリングを行った．
2.2 自由体積解析自由体積は既報 [4,9,10]の方法を用いて解析した．半径Rの球状プローブをゲスト分子として，高分子などのホスト分子に重なり合うことなく挿入することができる領域を，半径Rのプローブに対する自由体積(accessible volume; ACV)と定義する．ここでは，高分子を構成する各原子を，σ/2の半径を有する剛体球としてモデル化した．ここで，σはLennard-Jonesポテンシャルのサイズパラメータである．基本セルを0.2 Å間隔の格子点に分割し，それぞれの格子点に挿入することができる最大のプローブ半径Rmを求めた．Rm ≥ Rを満たす格子点が，半径Rのプローブに対するACVとなる．
また，ACVのクラスター解析を行うことにより，自由体積の結合性を調べた．ACV中の各格子点は，隣接格子と6つの結合を有する．隣接する2つの格子点が，ともに半径Rのプローブにアクセス可能である場合，それらは同一のクラスターに属するものとする．
2.3 3Dプリンターによる造形分子キャビティーの三次元モデルを作成するために，3D Systems社製Cube® 第二世代3Dプリンターを用いた．加熱された溶融樹脂材料がプリントチップから細く押し出され，プリントパッド上に層を形成する．各層が形成される度にプリントパッドが少しずつ下降し，新しい層が前の層の上に積み上げられていき，造形物の最後の層が一番上に形成されるまでこのプロセスが繰り返される．Cube® 3Dプリンターの特長を，以下に示す [11]．
・最大プリントサイズ 140 × 140 × 140 mm
・レイヤーの厚み 0.25 mm
・簡単に取り付け可能なプリントカートリッジ
・丈夫なABS及びPLA樹脂
・USB及びWi-Fi接続
立体モデル作成までの流れを，Figure 1に示す．上述の自由体積解析から得られた各格子点上のRmの値を構造格子データ(Fieldファイル)として出力し，可視化ソフトウェアMicroAVS [12]を用いて，Rm = R(Rの値については後述)となる面を"等数値面”として可視化した．これにより，半径Rのプローブに対するACVの表面が三次元表示される．得られた三次元形状データをSTL (Standard Triangulated Language)形式のファイルに出力した．必要に応じて，ソフトウェアnetfabb [13]を用いて，STLデータの編集，修正を行った後，3Dプリンター付属のCube® Softwareに読み込み，準備処理を行った．これにより，サポート構造(宙に浮いた部分やオーバーハング部を支える支柱)が生成され，スライス処理が行われる．処理したデータをCube形式で出力し，3Dプリンターに転送後，モデルを造形した．その後，造形物に付加された不要なサポート材料を，ラジオペンチやカッターナイフを用いて慎重に取り除いた．
Schematic diagram of modeling of molecular cavities by 3D printer.
分子キャビティー中に気体を挿入した場合については，一つの気体分子に注目し，初期位置からの変位のグラフよりジャンプする瞬間を特定し，軌跡ファイルからその瞬間の座標を抽出した．自由体積に関しては，前述と同様に，自由体積解析から得られたFieldファイルをMicroAVSで読み込み，STL形式で出力した．一方，気体分子に関しては，シミュレーションで得た原子の座標をPDBファイルとして出力し，これをソフトウェアUCSF Chimera [14]で読み込み，空間充填球モデルの三次元形状データをSTL形式で出力した．得られた二つのSTLファイルをnetfabbで同時に表示し，STLデータの編集，修正を行った．Cube®ソフトウェアで準備処理の後，データをCube形式で出力し，3Dプリンターでモデルを造形した．
樹脂に関しては，ABS樹脂は熱収縮による反りの発生が顕著であるため，PLA (ポリ乳酸)樹脂を使用した．また，形状によってはサポート構造が多く形成され，除去が困難となる．この場合，netfabbを用いて形状データを部分に分割し，それぞれをパーツとして造形した後，接着剤を用いて貼り合わせた．
なお，造形の詳細な手順については，電子付録「立体モデル造形の詳細な手順」に示した．
シミュレーションから得られた各結晶の単位胞の長さ，角度及び密度を，実験値と合わせてTable 2に示す．δe型とε型について，実験値と比較すると，誤差も少なく，構造はシミュレーションによって十分に再現されているといえる．S-I型は，b軸方向の応力280 MPa (∑yy = 16.0 kJ/molÅ)における値である．実験ではまだ報告されていないが，既報 [8]のシミュレーションとほぼ同じ構造が得られている．
| Form | a (Å) | b (Å) | c (Å) | γ (deg.) | d (g/cm3) |
| Simulation | |||||
| δe | 17.40 | 11.70 | 7.81 | 114.91 | 0.960 |
| ε | 16.16 | 21.42 | 7.96 | 90.00 | 1.004 |
| S-Ia) | 21.87 | 16.54 | 7.88 | 90.00 | 0.971 |
| Experiment | |||||
| δe[5] | 17.4 | 11.85 | 7.70 | 117 | 0.977 |
| ε[6] | 16.1 | 21.8 | 7.9 | 90.00 | 0.98 |
a) b軸方向の応力280 MPaにおける値
δe型について，シミュレーションから得られた結晶構造に対して，自由体積のクラスター解析を行った．この系では，ACVはベンゼン環程度のサイズの比較的大きなキャビティーと，これらをつなぐ細いチャネルからなっている．このチャネルより小さなプローブの場合はACVクラスターが系全体に連続的にひろがるが，プローブ半径がチャネル半径より大きくなると，ACVがチャネル部で分断され，クラスターが急激に小さくなる．このときの臨界的なプローブ半径Rcが，キャビティーをつなぐチャネルの半径に相当するといえる [4]．再現した結晶構造に関しては，Rc ≈0.60 Åであると判断できた．
プローブ半径RをRc付近としてACVのクラスター解析を行うことにより，気体のジャンプ拡散が生じる細いチャネル部の構造がより明瞭となる．Figure 2に，プローブ半径R = 0.6 ÅとしたときのACVクラスターの投影図を示す．これより，δe型結晶ではb軸方向及びac軸方向に，キャビティー間をつなぐ細いチャネルが存在していることがわかる．
Projections of one of continuous accessible clusters in the δe form calculated for a probe of radius 0.6 Å.
δe型と同様に，ε型及びS-I型についても，シミュレーションから得られた結晶構造に対して，自由体積の解析を行った．ε型とS-I型については，同じプローブ半径で構造を解析して比較するため，どちらもR = 0.7 Åとして解析を行った．このときのACVクラスターの投影図を，Figure 3に示す．Figure 3 (a)より，ε型ではc軸に沿って円筒型チャネルが存在している一方で，Figure 3 (b)より，S-I型ではc軸に沿ってジグザグ状のチャネルが存在していることがわかる．
Projections of one of continuous accessible clusters calculated for a probe of radius 0.7 Å.
δe型について，R = 0.6 Åのプローブを用いてACVを求め，立体モデルを作成した．b軸方向及びac軸方向にキャビティー間をつなぐように存在するチャネルの立体モデルをFigure 4に示す．紙面に水平な面で上下2つのパーツに分けて作成し，接着剤で貼り合わせたものである．
Channel models of the δe form crystal (left: b-axis direction, right: ac-axis direction).
Figure 4より，b軸方向よりもac軸方向につながるチャネルの方が太くて短いことがわかる．これまでの研究で，δe型中のキャビティー内の気体分子に関する優先的な拡散経路は，ac軸に沿って隣接するキャビティー間をつなぐチャネルを通る経路であることが報告されている [4,15]．Figure 4の結果は，キャビティー内の気体分子がb軸方向よりもac軸方向に拡散しやすいということと関連付けることができ，その差が十分に再現されていた．
ε型及びS-I型については，半径R = 0.7 Åのプローブを用いて自由体積解析を行い，立体モデルを作成した．ε型では，Figure 5 (a)のようにc軸に沿った円筒型チャネルのモデルを，S-I型では，Figure 5 (b)のようにc軸に沿ったジグザグ状のチャネルのモデルを再現することができた．基本的にFigure 3が立体モデルになっただけではあるが，実際に手に取って観察することにより，分子キャビティーの構造を直感的に理解できる意義は大きい．ε型のチャネルは円筒型であるが，所々b軸方向に飛び出た部分が見られる．応力を加えてS-I型に変形する際に，この部分が分かれてジグザグ状のチャネルを形成していることが考えられる．
Channel structure of s-PS crystals modeled by 3D printer.
S-I型結晶を用いることにより非常に効率よくCO2を分離できることが，著者(玉井)が行った分子シミュレーションにより予測されている [16]．N2やCH4と比べてCO2が非常に効率的に拡散することが示されており，その機構の解明が望まれる．そこで，ε型とS-I型について，分子キャビティー中に気体を挿入し，ジャンプ拡散の瞬間をとらえて三次元モデルを作成し，観察することにした．気体にはCO2及びN2を用いた．
ε型結晶中のある気体分子に注目し，初期位置からの変位を求めたところ，Figure 6のようになった．c軸方向の一次元拡散となり，また，比較的連続拡散に近いジャンプ拡散となっている．Figure 6の矢印で示した部分の瞬間的な結晶構造を抽出し，自由体積及び注目した気体分子の立体モデルを作成したところ，Figure 7のようになった．なお，用いた3Dプリンターでは2種類以上の樹脂を同時に用いることができないため，ACVと気体分子を別々に造形して，接着剤で貼り合わせている(詳細な手順に関しては，電子付録を参照)．
Displacement of gases in the s-PS ε form.
Instantaneous structure at the moment the gas molecule hops in the ε form, modeled by 3D printer.
Figure 6の気体の変位を考慮してFigure 7のモデルを観察すると，ε型のチャネルはc軸に沿った円筒型のものであるため，気体はチャネル内を上下に連続的に運動しており，その一瞬を再現できたことがわかる．また，気体分子が若干大きく表現されている(電子付録S2.1節参照)ことを考慮すれば，ε型のチャネルには気体が透過するのに十分な領域が存在していることが確認できる．
S-I型についても，ε型と同様の手順で気体を拡散させ，気体分子の変位を求めたところ，Figure 8のようになった．ε型と異なり，明確なジャンプ拡散となっていることがわかる．このジャンプの瞬間(図中の矢印で示した時刻)の構造を抽出し，自由体積及び注目した気体分子の立体モデルをそれぞれ造形し，接着剤で貼り合わせることにより，Figure 9に示すモデルを得た．
Displacement of gases in the s-PS S-I form.
Instantaneous structure at the moment the gas molecule jumps through a channel in the S-I form, modeled by 3D printer.
S-I型結晶では，気体分子程度の大きさのキャビティーがac軸方向の細いチャネルでジグザグに連結されている．Figure 8の気体の変位を考慮してFigure 9のモデルを観察すると，気体はチャネル内をac軸方向に跳躍的に運動しており(この間，b軸方向にはジャンプしない)，細いチャネル部分をジャンプする瞬間を再現できたことがわかる．また，気体分子が存在している位置では自由体積のモデルが非常に細くなっていることから，気体が透過するには，チャネルに対する気体分子の配向について，一定の条件が必要であることが考えられる．Figure 9のモデルにおいては，CO2，N2共に，気体分子の長軸が，チャネルがつながる方向と一致していることが確認できる．
本研究では，s-PS δe型，ε型及びS-I型結晶に対してMDシミュレーションを行い，得られた結晶構造に対して自由体積を解析し，3Dプリンターを用いてそれらの分子キャビティーの三次元構造モデルを作成した．
δe型に関しては，b軸方向及びac軸方向につながるチャネルのモデルを再現した．b軸方向よりもac軸方向につながるチャネルの方が太くて短いことがわかり，キャビティー内の気体分子がac軸方向に拡散しやすいということと関連付けることができた．ε型及びS-I型に関しては，それぞれのモデルを比較したところ，ε型で円筒型であったチャネルが，応力印加後のS-I型ではジグザグ状のチャネルに変化していることが確認できた．
また，ε型とS-I型については，分子キャビティー中に気体を挿入し，その拡散の様子を再現した．円筒型チャネルのε型では，チャネル内を振動する気体の様子を，ジグザグ状チャネルのS-I型では，ac軸方向につながる細いチャネル内を気体が跳躍する瞬間を再現できた．
今回，分子科学の分野で3Dプリンターを用いることで，分子キャビティーの複雑な構造をモデル化することができ，これまで把握することが困難であったその構造を，立体的・直観的に把握することができた．さらに，平面では比較できなかった，離れた位置にある構造や異なる構造を，モデル化することで並べて比較することができ，構造の解明につなげることができた．また，二酸化炭素や窒素などの気体がキャビティー中をジャンプ拡散する瞬間のモデルを作成し，高分子結晶中における気体の拡散機構に関して，自由体積の構造と気体の運動との関係性についての理解が深まった．今回のシミュレーションでは行わなかったδe型における気体の拡散についても，ε型，S-I型と同様な方法でモデル造形を行うことで，その解明につなげることができる．
課題としては，造形の精度の問題が挙げられる．本研究では，δe型のモデルに関して，平面で切断したものを接着剤で貼り合わせるという方法でモデルを作成したが，樹脂の熱収縮性のため，それぞれのパーツで若干の反りが生じ，接着面の端がきれいに接着しないという問題があった．δe型のモデルが熱収縮の影響を受けやすい形状であるためであると考えられるが，保温機能がついた3Dプリンターを用いることにより改善される可能性がある．また，気体拡散モデルにおいて，気体分子がプローブ半径に相当する分だけ大きく表現されている点についても，検討の余地がある．
これまで，分子シミュレーションで得た三次元構造は計算機のモニタ画面で観察していたが，これを手にとって観察できるようになったことは新鮮な体験であった．今後，タンパク質や低分子，他の高分子の系において，分子シミュレーションと3Dプリンターによる造形を融合することで，分子科学分野のさらなる発展に寄与することが期待される．
