原子间相互作用对无序大分子...和动力学影响的中子散射研究_韩泽华.pdf

原子间相互作用对无序大分子结构和动力学影响的中子散射研究韩泽华1,2,左太森1,2,马长利1,2,李雨晴1,2,3,程贺1,2*1.中国科学院高能物理研究所,中国散裂中子源,东莞 5238002.散裂中子源科学中心,东莞 5238003.中国科学院大学,北京 100049*通讯作者,E-mail:收稿日期:2022-12-30;接受日期:2023-02-07;网络版发表日期:2023-03-31国家自然科学基金(编号:22273046,U1932161)资助项目摘要不同单体在原子尺度的相互作用使得不同大分子体系在微观、介观的多尺度复杂结构及动力学行为有明显差异,从而进一步影响了体系的宏观性质.联用中子散射和计算机模拟,利用相同分子量和分子量分布的氘代大分子与氢化大分子具有相同分子结构、不同中子衬度,以及中子散射和计算机模拟宽广的时间、空间观察尺度,我们可以得到无序大分子体系最可几全原子结构,进而分析其从原子到纳米的多尺度空间结构,与从皮秒到微秒的多模式动力学行为的成因.近年来,我们使用该方法,从小分子到大分子的稀溶液、溶胀体系,从大分子熔体到玻璃态,成功分析了原子间相互作用对不同空间尺度结构和跨时间运动模式动力学行为的影响.本文介绍了这些典型的例子,希望将该方法推广到更广阔的研究领域,把大分子原子结构的多样性与多尺度的复杂结构和动力学有机地联系起来.关键词无序大分子中子全散射,准弹性中子散射,计算机模拟,全原子结构,多尺度结构和动力学1引言因为高分子材料的性能由组成其的单体成分和加工过程中形成的多尺度结构共同决定,所以高分子领域的物理学家和化学家分别采用独特的方法来提高材料的性能.物理学家的解决方式是建立具有一般性的理论来描述体系中结构、动力学行为和宏观性质的关系,并将其应用到未知体系中预测性能和表现.比如,Flory1的平均场理论就建立在粗粒化的格子模型上,通过引入不可压缩假定,描述了高分子的基本相分离过程;de Gennes的标度理论同样建立在粗粒化的串滴(blob)模型上,描绘了温度blob和浓度blob的标度变化等2.而高分子化学家的解决方案则是合成由不同单体组成的新型高分子,利用不同单体间复杂的相互作用组装成各种各样的多尺度结构.在两种方式相互联系的诸多桥梁中,经验参数和经验公式发挥了重要的作用.例如,依据溶度参数预测高分子在不同溶剂中的溶解能力,从而选择不同溶引用格式:Han Z,Zuo T,Ma C,Li Y,Cheng H.The influence of interatomic interaction on the structure and dynamics of disordered macromolecules:a neutronscattering study.Sci Sin Chim,2023,53:678692,doi:10.1360/SSC-2022-0257 2023 中国科学杂志社中国科学:化学2023 年第 53 卷第 4 期:678 692SCIENTIA SINICA C聚合物结构与性能专刊专题论述剂使大分子自组装3,4;依据侧基的空间位阻来推测高分子膜的气体选择能力5;依据Vogel-Fulcher-Tam-mann(VFT)方程来计算高分子玻璃化转变过程的动力学特征弛豫时间等6.无法定量这些经验参数与经验公式的主要原因是研究手段的匮乏:传统的研究手段无法在高分子材料本体中实时观测其全原子结构及动力学,只能通过微观、介观、宏观尺度上的一些结构与性能或者光谱的变化推测多尺度结构成因的原子图像7,8,由于这种间接观察的方式不直观,所以研究者很难从单独个体特性出发,总结出系统性的结论.中子散射与计算机模拟相结合,可以有效地解决这个问题.一方面,二者都可以覆盖从埃到微米、从皮秒到微秒的多尺度范围;另一方面,二者又可以互补不足,互相支撑.中子的优点主要表现在:强穿透性,可以在复杂的样品环境下实现对高分子体系结构(变化)的实时观测;氘代技术,可以选择性观察复杂体系中不同组分的结构和动力学等9,10.中子散射的主要缺点是不直观:作为一种倒空间的系综平均,中子散射无法给出体系中微观过程的直观物理图像;而计算机模拟最大的优势就是直观,通过模拟结果可以直接确定体系实空间的结构与运动:通过使用全原子或者粗粒化力场对实验体系进行多尺度模拟,直接与中子散射实验结果对比,分析多尺度结构和动力学的成因,进而拓展计算,从而得到许多无法直接从中子散射实验中得到的细节,解释经验参数背后的物理.迄今为止,该研究方法已经取得众多成果,如英国散裂中子源(ISIS)的Soper等1116使用无序小分子中子全散射结合蒙特卡罗模拟的方法分析了多种小分子溶液的结构;美国麻省理工的Chen等17使用准弹性中子散射与计算机模拟结合的方式研究分析了水在220 K附近时的动力学转变过程;西班牙材料物理中心的Colmenero等18使用小角中子散射、准弹性中子散射结合分子动力学模拟的方式系统研究了线性共聚高分子中的特征结构、动力学行为与单体间关系;加拿大不列颠哥伦比亚大学的Cullis等19利用小角中子散射与计算机模拟结合的方式得到了包裹siRNA的脂质体结构;瑞典隆德大学的Wacklin等20利用分子动力学模拟分析中子反射的实验结果表征了生物脂质体膜的特征结构;埃克森美孚石油的Lpez-Barrn等21使用小角中子散射、广角X射线衍射及分子动力学模拟结合的方式分析了聚乙烯基联苯(poly(vinylbiphenyl),PVBP)中的局部特征结构并得到了实空间图像;上海交通大学的洪亮等22,23使用准弹性中子散射与计算机模拟结合的方式分析了低温下水的行为与蛋白质行为的关系:辅仁大学的Yang等24使用小角中子散射结合计算机模拟的方式系统研究了多种蛋白质在溶液中的结构等.本课题组结合无序大分子中子全散射、准弹性中子散射与计算机模拟,开发了一套解析无序大分子最可几全原子结构和动力学的方法.从Soper教授20世纪80年代开创的实验势结构精修(empirical potentialstructure refinement,EPSR)模拟无序液体全散射数据的成功经验出发16,引入GPU加速和并行计算25,从多组分小分子混合体系开始,随后在大分子的稀溶液、浓溶液、溶胀以及熔体和玻璃体系中进行了实验验证,得到了大分子的最可几全原子结构与动力学图像,分析了大分子从最可几全原子位置到整体构象的多尺度特征结构,以及从皮秒到微秒的多模式(mode)运动行为2632.该方法可以有效解释不同单体、不同原子间相互作用对多尺度空间结构和动力学的影响.本文先介绍了所用研究方法的基本原理,然后介绍了使用该方法在小分子溶液、大分子溶液、大分子熔融态和玻璃态体系中对结构与动力学行为研究得到的成果.2理论基础在中子全散射中,微分散射截面F(q),即qdd()代表单个散射体将中子散射到立体角d内的概率33,其与散射体内粒子的空间结构有如下关系34:()()c bc c b brgrqrqrrqdd()=+24()1sind(1)_incoh20+2其中q为散射矢量,、代表不同的原子类型,c代表该原子数目占体系原子总数的比值,b代表散射长度,为体系中全部原子的数密度,g(r)为径向分布函数,代表以粒子为中心时,粒子数密度对平均密度的比值沿距离r的分布情况35,计算方式为:中国科学:化学2023 年第 53 卷第 4 期679grrrrr()=1()4(2)2其中代表原子的数密度.当使用氘代技术改变体系中氢(H)与氘(D)的比例后,同一结构可得到多条不同的散射曲线,其对应的实空间结构完全相同25,即:FFwwwwHHqq().()=.(3)MNMMNN111111其中F(q)为不同氢氘比下体系的散射结果,w为粒子间相互作用的散射权重因子:()wc c b b=2(4)ijijijijH为粒子间相互作用:()Hrgrqrqrrq()=4()1sind(5)ijij0+2i与j由w矩阵与H矩阵中对应元素的列号与行号决定.根据式(1),如果各向同性的无序大分子体系中有N种不同的原子,就存在N(N+1)/2即(CN+N2)种不同的g(r),即如果想解出所有的g(r),理论上需要有N(N+1)/2个不同衬度的独立实验.在实际研究中,我们可以通过调整样品中不同组分的氘代情况,得到尽可能多的独立散射结果,同时利用计算机模拟预测补充实验数据的不足:使用同一体系不同衬度下的多个I(q)数据与对应组分的计算机模拟结果进行对比,得到体系中所有原子的最可几位置,这可以使模拟结果的散射数据验证过程更加可靠.对于小分子体系的模拟,我们主要在Soper教授的EPSR方法的基础上发展了EPSR+的程序,添加了GPU加速与并行计算来进行拟蒙特卡罗(逆蒙卡)模拟25.逆蒙卡模拟在传统蒙特卡罗(蒙卡)模拟的基础上,以散射结果为目标,通过迭代优化模拟体系中使用的势能函数,从而优化计算机模拟最终得到的结构.该模拟方法与小分子体系的中子全散射实验方法的匹配度非常好,然而在大分子体系中,由于大分子特征尺寸较大,其结构信息是通过小角散射区域内的数据体现,往往没有明显的特征峰,因此小角散射的实验数据较难通过傅里叶变换得到准确的势能函数结果,这也导致迭代过程比较不稳定.在大分子体系中使用较多的是分子动力学(mole-cular dynamics,MD)模拟,即通过设置粒子间的相互作用关系,计算其受力情况及相应的运动状态来进行体系的模拟.粒子间的相互作用参数,即力场的设置,对于模拟结果具有决定性的作用.一般来说,在分子动力学模拟中会涉及的粒子间相互作用包括成键、键角、二面角以及非键相互作用36.较常用的成键相互作用能与键角相互作用能的计算方式均为简谐式,即:EK rr=()(6)rbondeq2EK=()(7)angleeq2其中脚标eq表示平衡位置,K为作用参数,二面角相互作用的计算略有不同,为EVVV=21+cos()+21cos(2)+21+cos(3)(8)torsion123其中Vi为对应项的傅里叶系数,为二面角.粒子间的非键相互作用一般包含两项,即范德华作用力与静电力,其势能函数为:Efqq errr=+4(9)n bijijijijijijijij.2121266其中f为权重因子,q代表粒子的电荷数,e为电子电荷,r为粒子间距离,与为范德华相互作用参数.通过式(6)(9)可以计算体系中所有粒子间的相互作用,进一步可以得到粒子的受力情况:F rU rr()=d()d(10)再利用牛顿力学原理计算粒子的运动状态.粒子的初始速度可以自行设置合理数值,也可以通过模拟体系所处的温度根据泊松分布进行分配.我们可以根据已有的可靠力场进行模拟,得到体系的最可几全原子结构,然后根据式(1)进行计算,并与中子散射实验结果进行比较:如果实验结果与同一体系多种不同衬度下样品的实验结果都吻合得很好,则可以认为得到了体系的最可几全原子结构.粗粒化模拟的原则和全原子模拟基本相似,需要预先设定体系的力场参数再进行模拟,力场所包含的势能相互作用基本上与式(6)(9)相同.为了能保证粗粒化模拟的结构与中子散射结果相吻合,可以利用迭代玻尔兹曼反演方法,通过全原子模拟的结果得到粗韩泽华等:原子间相互作用对无序大分子结构和动力学影响的中子散射研究680粒化体系中的力场参数,从而保证粗粒化体系在单体尺度上的结构与全原子的结构完全相同37.该算法的基本思路为通过粗粒化结果与全原子结果间的差异迭代势能函数,最终得到更精确的势能函数:U rkTgr()=ln()(11)0targetUrU rkTgrgr()=()+ln()()(12)ii+1calculatedtarget其中U为势能函数,脚标i与i+1代表迭代中的轮次,g(r)在这里代表的是两体相互作用所涉及的分布函数,既可以是键长、键角分布,也可以是非键相互作用粒子间的相对位置分布,target为全原子模拟的结果,calcu-lated为粗粒化模拟计算得到的结果.当然,MD模拟也有自己的局限性,其中