|
2017年05月02日 19:58:51
看到一篇科技论文的报道纳米粒子反常受限扩散研究获得进展, 是有关扩散的, 就顺便整理一下与扩散有关的理论知识, 这些知识在利用分子动力学模拟研究扩散现象时可能会用得上.
扩散模式的不同分类使用MD计算粒子的扩散系数时, 一般都是基于粒子进行布朗运动, 遵循简单的扩散模型. 实际上, 根据外界限制条件的不同, 粒子的扩散有多种模式, 不同扩散模式下粒子的均方位移MSD(mean square displacement)与时间t的关系很不一样. 根据MSD的表现, 我们可以推测出扩散的不同机理. 研究粒子, 如量子点在细胞中的运动时, 人们总结了粒子几种可能的扩散模式. 下面是论文Confined Lateral Diffusion Of Membrane Receptors As Studied By Single Particle Tracking (Nanovid Microscopy). Effects Of Calcium-induced Differentiation In Cultured Epithelial Cells中给出的总结.
静态模式 Stationary mode
粒子基本不运动, 扩散系数接近零.
简单(或布朗)扩散模式 Simple diffusion mode
粒子做简单的布朗扩散运动
MSD(Δt)=nDΔt, <span class="MathJax" id="MathJax-Element-2-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">n" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">n 为粒子运动空间的维数
定向扩散(或传输)模式 Directed diffusion mode (transport mode)
粒子在做随机扩散的同时沿某一方向以恒定的速度漂移
<span class="MathJax" id="MathJax-Element-3-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">MSD<mo stretchy="false">(<mi mathvariant="normal">Δt<mo stretchy="false">)=nD<mi mathvariant="normal">Δt+v2<mi mathvariant="normal">Δt2" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">MSD(Δt)=nDΔt+v2Δt2
限制扩散模式 Restricted diffusion mode
粒子在做布朗扩散时, 其运动空间被限制在某一范围内, <span class="MathJax" id="MathJax-Element-4-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">0≤x≤Lx" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">0≤x≤Lx , <span class="MathJax" id="MathJax-Element-5-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">0≤y≤Ly" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">0≤y≤Ly , <span class="MathJax" id="MathJax-Element-6-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">0≤z≤Lz" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">0≤z≤Lz . 粒子的运动等价于在无限深方势阱中的布朗扩散. 细分起来又可分为两类: 受限模式confinement model和系索模式tethering model. 但仅根据轨迹无法区分这两种模式.
<span class="MathJax" id="MathJax-Element-7-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML"><mtable columnalign="right left right left right left right left right left right left" rowspacing="3pt" columnspacing="0em 2em 0em 2em 0em 2em 0em 2em 0em 2em 0em" displaystyle="true">⟨x2⟩<mo stretchy="false">(t<mo stretchy="false">)=<mrow class="MJX-TeXAtom-ORD">Lx26−<mrow class="MJX-TeXAtom-ORD">16Lx2π4<mo movablelimits="false">∑<mrow class="MJX-TeXAtom-ORD">n=1<mo stretchy="false">(odd<mo stretchy="false">)<mi mathvariant="normal">∞<mrow class="MJX-TeXAtom-ORD">1n4exp⁡[−<mrow class="MJX-TeXAtom-ORD">12<mo stretchy="false">(<mrow class="MJX-TeXAtom-ORD">nπσxLx<mo stretchy="false">)2t],σx2=2Dx⟨y2⟩<mo stretchy="false">(t<mo stretchy="false">)=<mrow class="MJX-TeXAtom-ORD">Ly26−<mrow class="MJX-TeXAtom-ORD">16Ly2π4<mo movablelimits="false">∑<mrow class="MJX-TeXAtom-ORD">n=1<mo stretchy="false">(odd<mo stretchy="false">)<mi mathvariant="normal">∞<mrow class="MJX-TeXAtom-ORD">1n4exp⁡[−<mrow class="MJX-TeXAtom-ORD">12<mo stretchy="false">(<mrow class="MJX-TeXAtom-ORD">nπσyLy<mo stretchy="false">)2t],σy2=2Dy⟨z2⟩<mo stretchy="false">(t<mo stretchy="false">)=<mrow class="MJX-TeXAtom-ORD">Lz26−<mrow class="MJX-TeXAtom-ORD">16Lz2π4<mo movablelimits="false">∑<mrow class="MJX-TeXAtom-ORD">n=1<mo stretchy="false">(odd<mo stretchy="false">)<mi mathvariant="normal">∞<mrow class="MJX-TeXAtom-ORD">1n4exp⁡[−<mrow class="MJX-TeXAtom-ORD">12<mo stretchy="false">(<mrow class="MJX-TeXAtom-ORD">nπσzLz<mo stretchy="false">)2t],σz2=2Dz6D=2Dx+2Dy+2Dz" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">⟨x2⟩(t)=L2x6−16L2xπ4∑n=1(odd)∞1n4exp[−12(nπσxLx)2t],σ2x=2Dx⟨y2⟩(t)=L2y6−16L2yπ4∑n=1(odd)∞1n4exp[−12(nπσyLy)2t],σ2y=2Dy⟨z2⟩(t)=L2z6−16L2zπ4∑n=1(odd)∞1n4exp[−12(nπσzLz)2t],σ2z=2Dz6D=2Dx+2Dy+2Dz
受障模式 Obstacle-impeded diffusion mode
粒子进行自由扩散时受障碍物限制, 障碍物可以是固定的, 或者有一定的移动能力. 在这种情况下, 长程扩散会减弱但仍大于零. 这种模式比较难与限制模式区分开来. 对各向同性的碰撞几率,
<span class="MathJax" id="MathJax-Element-8-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">MSD<mo stretchy="false">(<mi mathvariant="normal">Δt<mo stretchy="false">)=4D<mi mathvariant="normal">Δt+A+Bln⁡<mo stretchy="false">(C<mi mathvariant="normal">Δt<mo stretchy="false">)" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">MSD(Δt)=4DΔt+A+Bln(CΔt)
大多数情况下MD研究的是简单扩散模式, 在此假定下求扩散系数的步骤是先算MSD, 然后拟合MSD线性部分的斜率, 再除以6(二维体系除以4)就是扩散系数. 但在拟合的时候, 有一个选取拟合数据范围的问题, 也就是确定用哪段时间范围内的MS进行拟合. 很显然, 使用不同时间段的数据, 拟合结果会有所不同. 拟合时不能使用MSD数据的起始部分, 因为起始部分属于扩散弛豫过程, 除非专门研究这个过程, 拟合的初始时间要大于扩散弛豫时间; MSD数据的最后部分也不能用, 因为计算MSD时, 关联时间越大, 数据点越少, MSD的误差越大, 所以拟合时只能取数据中间接近线性的一部分.
如果是普通扩散过程, 在开始的一小段时间内MSD是关联时间的二次函数, 代表无障碍的定向扩散. 随着关联时间的增加, MSD会很快过渡到一次函数阶段, 代表正常扩散. 这个一次函数区域一般来说就是计算扩散系数的最佳区域.
由于关联时间越大, 涨落越大, 所得MSD的误差也越大, 所以拟合的最大关联时间也不是越大越好. 此外, MSD曲线的光滑程度与模拟时间成正比, 一般模拟时间会取最大关联时间的10倍或更大. 也可以多次模拟求平均值.
具体怎么决定拟合的关联时间范围呢? 一种方法是选取不同的最大关联时间, 如1 ps, 5 ps, 10 ps, 50 ps, 100 ps等作出MSD随最大关联时间变化的曲线, 从而确定体系的特征扩散弛豫时间. 更好的方法是计算动扩散系数(RDC, running diffusion constant), 根据它的变化来确定关联时间的拟合范围
动扩散系数曲线要趋近于水平线才算收敛. 对于普通扩散, 随关联时间的增加, 动扩散系数会从零增加到一个定值. 对于气体, 增加过程一般是单调的; 对于稠密液体, 可能会有局部涨落. 当然,单次模拟数据会有较大涨落, 我们很难通过一次模拟结果就判断出这个定值是多少, 所以需要多次模拟取平均.
具体做法就是, 使用 MSD/6t 对t作图, 这样数据成正比的线性部分应该是一条上下稍有波动的水平线, 很容易看出来. 我们可以简单地将这段时间内扩散系数的平均值模拟的扩散系数, 也可以对此时间段内的数据进行拟合求出扩散系数. 这两种方法得到的扩散系数值应该差不多, 但建议采用后一种方法, 虽然麻烦点, 但标准
这种方法的一种变形是做双对数 log(MSD)-log(t) 图, 选取其斜率尽可能接近1的一段求扩散系数. 这很容易理解
缺点在于单凭肉眼很难判断数据斜率为1的部分.
更复杂的处理方法, 可以采用随机抽样一致性算法(RANSAC)或者非线性拟合方法来自动确定最佳拟合值, 其大致思想是在指定的条件下, 寻找一条直线, 使得距离直线一定范围内的点数最多. 但这些方法使用麻烦, 所以文献中较少看到.
MD计算设置系综选用: 大多用NVT, 但NPT也可以, 不过要注意NPT时盒子标度导致的位置变化
MSD计算属于相关函数计算, 需要大量的数据点, 因此, 需要长的模拟时间以及短的输出间隔. 只有模拟时间要足够长才能保证扩散系数收敛.
扩散系数
单位 [D] = nm2/ps = 10-6 m2/s = 10-2 cm2/s = 103 10-5 cm2/s = 103 10-9 m2/s
1 cm2/s = 100 nm2/ps
298.2 K时, 水的扩散系数 Dwat=(2.30 ± 1.5%) x 10-9 m2/s=(2.30 ± 1.5%) x 10-5 cm2/s.
其他一些分子的扩散系数, 可查阅CRC物理化学手册.
gmx msd计算时可以设定拟合时间, 计算后会输出相关信息以及扩散系数, 但如果你不能确定拟合时间范围的话, 给出的数值最好不要直接使用.
理论上讲, 当时间足够长时, 粒子的运动可以近似看作随机游走, 那么MSD与关联时间成线性关系. 在三维情况下, <span class="MathJax" id="MathJax-Element-11-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">MSD<mo stretchy="false">(t<mo stretchy="false">)=6Dt" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">MSD(t)=6Dt . 对上式两边取对数作图, 斜率为1. 如果计算MSD的时间不够长, 你会发现MSD与时间呈非线性关系 <span class="MathJax" id="MathJax-Element-12-Frame" tabindex="0" data-mathml="<math xmlns="http://www.w3.org/1998/Math/MathML">MSD=Dαtα<mo stretchy="false">(a<1<mo stretchy="false">)" role="presentation" style="margin:0px;padding:0px;display:inline;line-height:normal;text-align:left;word-spacing:normal;word-wrap:normal;float:none;direction:ltr;max-width:none;max-height:none;min-width:0px;min-height:0px;border:0px;position:relative;">MSD=Dαtα(a<1) . 这时对两边取对数作图斜率不为1. 这种情况称为sub-diffusive行为. 直观的解释是粒子的游走可以近似看作是在其他粒子形成的”溶剂笼”里活动, 粒子长时间运动的MSD可近似看作粒子在各个溶剂笼间的活动, 因此与时间呈线性关系, 而粒子短时间运动的MSD可看作是粒子还没有走出第一个溶剂笼, 因此是呈非线性的sub-diffusive行为. 为了克服sub-diffusive行为对计算扩散系数带来的影响, 通常是计算足够长时间的MSD, 将MSD/t或d(MSD)/dt对时间作图, 考察MSD斜率随时间的变化, 当斜率达到稳定值时就可以用来计算扩散系数了.
相关文献
Diffusion And Subdiffusion Of Interacting Particles On Comblike Structures
O. Bénichou, P. Illien, G. Oshanin, A. Sarracino, R. Voituriez; Phys. Rev. Lett. 115(22):, 2015; 10.1103/physrevlett.115.220601
Discriminating Between Anomalous Diffusion And Transient Behavior In Microheterogeneous Environments
Alexander M. Berezhkovskii, Leonardo Dagdug, Sergey M. Bezrukov; Biophysical Journal 106(2):L09-L11, 2014; 10.1016/j.bpj.2013.12.013
樊哲勇, 用分子动力学模拟计算扩散系数
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-9-26 06:20
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社