熔融Cu57团簇在急冷过程中弛豫和局域结构转变的分子动力学研究*

樊沁娜 李蔚 张林

(东北大学理学院，沈阳110004)

(2009年5月11日收到;2009年7月3日收到修改稿)

熔融Cu57团簇在急冷过程中弛豫和局域结构转变的分子动力学研究*

樊沁娜 李蔚 张林†

(东北大学理学院，沈阳110004)

(2009年5月11日收到;2009年7月3日收到修改稿)

采用基于嵌入原子方法的正则系综分子动力学研究熔融Cu57团簇在急冷过程中的弛豫及其局域结构变化.通过对弛豫过程中均方位移、非相干中间散射函数和非Gauss参数三种函数和原子键对随急冷温度不同所发生变化的分析表明，在经过短时间的原子剧烈运动后，急冷温度极大地影响着团簇内原子结构弛豫过程.急冷温度较高时，原子在经历短时间剧烈运动的β弛豫后，进入α弛豫区后以扩散运动为主，随后原子运动表现为非扩散性的原子局域结构重排，团簇内没有出现明显的成核结构.随着温度的降低，原子局域结构的变化在经过短时间原子剧烈运动的β弛豫后，在α弛豫区原子运动表现为扩散性运动，并出现一定数量的不稳定二十面体结构.当急冷温度很低时，在进入α弛豫区后，团簇结构变化逐渐表现为非扩散性原子局域结构重排，形成相当数量的稳定成核二十面体结构.

团簇，分子动力学，计算机模拟，表面

PACC:3640B，7115Q，6185，6800

1. 引言

货币金属金、银和铜都具有d10s1价电子结构，近些年来由这三种金属所组成的团簇，因其具有迥异于块体材料的丰富多样的物理与化学性质日益受到研究者们的广泛关注，并在生物及生物物理上具有广泛的应用［1—5］.在通过物理方法实验制备金属团簇的过程中，大多采用气相或离化的金属在不同方式下快速冷却凝结得到，尽管实验观测到结构各异的团簇结构，但对于包含几十个原子的小尺寸团簇在凝结过程中原子是如何运动及其所导致的结构变化目前尚知之甚少.了解这些小尺寸团簇在冷却凝结过程中的弛豫过程、成核和局域结构的变化，对于预测和控制团簇尺寸进而实验合成具有特定结构及性质的纳米颗粒具有指导意义.

为了研究金属熔融液体在过冷时所发生的结晶转变，1984年由Bengtzelius［6］等以及Leutheusser［7］在有关液态理论的基础上提出了MCT(mode-coupling theory)理论，这种理论认为结晶转变过程是一个动力学转变，在发生凝结的液体中存在两种弛豫过程，慢的α弛豫和快的β弛豫，α弛豫对应结构弛豫过程，该过程受温度影响较大，β弛豫是二级弛豫过程，该过程受温度影响较小.通过研究一些随时间变化相关函数的变化，MCT理论提出当冷却温度接近某一温度时，会在β弛豫末期也就是α弛豫前期出现一个过渡区间.应用该理论，研究者们已经对普通液体、金属熔融液体和较大的金属团簇在冷却过程中的弛豫特性作了较多的探讨［8—14］.但是对于包含几十个原子的小尺度自由表面金属团簇，相关的研究则尚未见到.这些小团簇由于它们具有很高的比表面(即表面积与体积之比)，使得它们随温度改变所发生的结构转变过程与大尺寸团簇和块体相比，展现出不同的特点.其中作为一种重要货币金属团簇的铜团簇以其在纳米催化剂、数据存储介质、医药诊断以及生物传感器等方面的潜在应用价值，备受研究者们的关注［15—20］.

本文采用基于嵌入原子方法(embedded atom method，EAM)的正则系综分子动力学(molecular dynamics，MD)，通过均方位移(mean squaredisplacement，MSD)、非相干中间散射函数(incoherent intermediate scattering function)和非Gauss参数(non-Gaussian parameter，NGP)等三种随时间变化的函数研究熔融Cu57团簇在被急冷到不同温度下的弛豫行为，通过键对分析(pair analysis，PA)方法给出了在不同急冷温度下团簇局域结构的变化.

2. 模拟方法

原子间的相互作用势采用Mei等［21］提出的EAM形式，应用该形式所给出势的分子动力学计算已用于模拟包含几十个原子小尺寸铜团簇在温度变化过程中的结构转变［15，17，19，20］.在EAM形式中，体系的总能量Et为

其中，ρi是其他原子在原子i处产生的电荷密度，F(ρi)是将原子i嵌入到电子密度为ρi处所具有的能量，rij是原子i和原子j之间的距离，(rij)是原子i和原子j之间的两体势.模拟采用正则系综(NVT)分子动力学方法，计算中时间步长取为1.6×10－15s.计算的几个函数如下，

其中〈·〉表示对于统计时间步的平均值.对分布函数g(r)给出了在相同原子密度下原子随机分布时，找到相距为r的一对原子的概率.V是所模拟体系的体积，N是体系中的原子数.〈r2(t)〉是均方位移函数MSD，其中ri(t)和ri(0)是第i个原子分别在时刻t和0时的位移.Fs(q，t)表示非相干中间散射函数，是von Hove函数自相关部分的Fourier变换，式中的波向量q选取为液态对分布函数的第一峰.α2(t)是非Gauss参数NGP，〈r2(t)〉是均方位移，〈r4(t)〉是平均四次方位移.非高斯参数被用于描述过冷液体的动力学非均匀特性，当液体是均匀分布并处于平衡态时，α2(t)接近于零.这里均方位移函数MSD，Fs(q，t)和α2(t)都是随时间变化描述团簇弛豫特性的函数.

研究中采用包含四个整数的键对分析技术研究急冷过程中原子局域结构的变化［22，23］.第一个整数表示一对原子是否为设定截断距离内的近邻原子(或成键原子对)，当它为1时表示成键.第二个整数表示某一成键原子对的共有近邻原子数，第三个整数表示这些共有近邻原子间的成键数目.第四个为附加数，表示某一特定的原子键对.例如，1551对应两个成键原子共有五个近邻原子，且这五个近邻原子形成具有五个键的五边形，该键对是二十面体结构的特征键对.1421和1422键对分别用于表征面心立方和密排六方的局域结构.

模拟中首先构造一个20a0×20a0×20a0的包含32000个原子的Cu晶体(a0=0.3615nm为Cu的晶格常数).20a0×20a0×20a0为模拟中分子动力学元胞的体积，在这个Cu晶体块体中截取出原子数目为57个原子的孤立自由表面铜团簇，将这个团簇直接升温到1800K，得到处于熔融状态的团簇.图1显示了这个团簇在1800K时的对分布函数曲线.由图中可见，1800K时Cu57团簇的对分布函数曲线表现出典型的液态结构特征，即对分布函数仅明显存在两个已经展宽的峰.在图1中右上侧内插图显示了该团簇在1800K时处于熔融无序状态的原子结构图.

图1 Cu57团簇在1800K时的对分布函数内插图为原子结构图

3. 结果分析

在获得处于熔融状态的Cu57团簇后，将其从1800K急冷到1200，1100，1000，900，800，700，600，500，400，300和200K时，图2给出了该熔融团簇被急冷到不同温度时的原子平均能量.由图中可见，Cu57团簇能量的降低可以分为三个阶段:在急冷温度高于1000K时，能量以较大的斜率近直线下降;当急冷温度低于500K后，能量以低于高急冷温度区间斜率的近直线形式降低;在500到1000K的急冷温度区间内，能量的变化呈现为振荡下降的形式.

图2 Cu57熔融团簇急冷到不同温度的原子平均能量

图3(a)，(b)和(c)分别显示了熔融铜团簇在被急冷到1200，600和200K时均方位移(MSD)、非高斯参数(NGP)和非相干中间散射函数Fs(q，t)随时间的变化.由图3(a)可见，当急冷温度为1200K时，MSD曲线表现为典型的液态行为.即在短时间0.1 ps的β弛豫区内，由于原子的集体剧烈运动，MSD以〈r2〉∝t2的形式迅速增加.随着保温时间的增加，原子的运动表现为扩散行为，进入α弛豫区，这时〈r2〉∝t.当时间增加到50 ps后，MSD的增速进一步放缓，这表明此时团簇内结构变化以非扩散性原子局域结构重排为主.这里，扩散性原子运动是指原子挣脱了其周围原子的束缚，此时该原子与周围原子构成的局域结构也发生了改变，非扩散性原子局域重排是指虽然某一原子运动，但它仍处于其周围原子的束缚中，即该原子与其周围原子构成的局域结构没有发生改变，从而在整体上看来是团簇内这些原子局域结构之间位置的重新排布.随着急冷温度的降低，MSD曲线出现了不同于较高温度时的形式.当急冷温度为600K时，在短时间0.1 ps内的β弛豫区内，MSD与时间的关系仍表现为〈r2〉∝t2，随后增速变缓，在经过这个缓慢变化区到0.2 ps后进入α弛豫区，MSD继续增加并与时间具有〈r2〉∝t的关系.随着急冷温度进一步降低到200K，进入MSD增速变缓区的时间缩短到0.08 ps，并在经过一个曲线缓慢增加的阶段后，在50 ps后进入一个平台区.如图3(a)200K时的MSD曲线所示，这个平台的高度要远低于高温时的MSD曲线.这说明当急冷温度很低时，每个原子的热运动降低，其移动范围变小，在团簇内原子间出现了较多的局域有序结构.另外根据Einstein关系，〈r2〉=6Dt，即当时间足够长时，自扩散系数D可以由MSD曲线的斜率计算得到，很明显，MSD出现平台意味着由于原子局域结构重排使得D趋于零.这样根据MSD曲线，我们就能确定在什么时间上过冷熔融团簇内开始了成核的转变过程.由图中可见，当急冷温度为600K时，过冷团簇内没有发生稳定的成核或者说团簇内局域有序结构还不能够较长时间的存在，而当急冷温度降低到200K时，过冷团簇在50 ps后出现了稳定的成核，这时存在了较多的局域有序结构.这种转变的进一步证据可由非相干中间散射函数Fs(q，t)和非高斯参数(NGP)显示出来.

图3 熔融Cu57团簇被急冷到1200，600和200K时的均方位移、非相干中间散射函数和非Gauss参数随时间的变化(a)均方位移MSD，(b)非相干中间散射函数Fs(q，t)，(c)非Gauss参数α2(t)

由图3(b)的非相干中间散射函数Fs(q，t)可见，Fs(q，t)曲线从Fs(q，0)=1开始，在短时间内即开始下降，随着温度的不同，该函数下降的行为也表现出差异.在急冷到较高温度1200K时，Fs(q，t)很快衰减为零，这表明此时团簇内的原子运动仍遵循通常液态的行为.随着急冷温度降低到600K时，团簇中的原子会出现结构重组，这时与MSD曲线的变化相对应，在经过一个快速下降的过程后，Fs(q，t)的变化出现一个较缓的下降过程，随后再次下降趋于零，这时非相干中间散射函数趋于零的时间明显比较高急冷温度时要长.当急冷温度下降到200K时，Fs(q，t)曲线在0.1 ps时即进入曲线以较缓坡度下降区间，并在50 ps后出现了一个平台，这说明随着急冷温度的降低，团簇原子已不会以液态原子的形式运动，并发生结构重排，这时Fs(q，t)不会趋于零而是出现一个平台.

类似的情形也出现在随时间变化的非Gauss参数变化上.如图3(c)所示，在急冷温度为1200K时，NGP没有表现出明显的峰.但当急冷温度降低到600和200K时，NGP曲线的各峰出现了升高，这表明被冷却系统结构变化的复杂程度增加，并且当急冷温度为200K时，NGP曲线在温度下降过程中出现了非零的平台，平台的出现说明此时团簇内原子已可以长时间处于局域有序的平衡位置.这里需要指出的是，这三个函数在急冷到不同温度时随时间的变化基本符合，即当急冷温度较高时，团簇内原子所进行的β弛豫和α弛豫如同通常的简单液体.随着急冷温度的降低，除了通常的β弛豫外，在α弛豫区内会出现一个短时间的非扩散性原子运动行为，但随后原子们的运动又以扩散性的原子运动为主.当急冷温度很低时，团簇内原子在经过集体剧烈运动的β弛豫后，α弛豫区内逐渐表现为非扩散性的原子局域结构重排.这样上述三个函数的变化表明所研究团簇在急冷过程中会出现原子集体运动，即当一个被束缚的原子要移动加快或减慢时，它的近邻原子也会如此，这必然导致所研究体系局域结构的变化.

图4 急冷温度为1200，600和200K时原子键对相对数目随时间的变化右侧图为其所对应的随时间变化的原子结构图. (a)1200K，(b)600K，(c)200K

图4(a)，(b)和(c)分别给出了团簇内原子在被急冷到1200，600和200K温度时的1421，1422和1551三种原子键对相对数目随时间的演变，这里的原子键对相对数目由某一键对数目除以键对总数.图的右侧为随时间变化的原子结构图.在模拟中的不同温度处，还会出现其他表征菱形对称结构的1201，1301和1311等键对以及表征液态和非晶态结构的1541和1431键对等.图4之所以只选取1421，1422和1551三种键对是考虑到在实验和模拟中发现大部分的铜团簇具有面心立方结构或二十面体结构，特别是包含五十几个原子的小尺寸团簇具有基于二十面体的几何形状.由7个原子组成的二十面体基元具有1551键对特征，随着由更多的原子组成稍大的二十面体，会出现较多数量的1422键对，并且当团簇尺寸更大时，所形成的二十面体还会有相当数量的1421键对.由图4中可见，当急冷温度为1200K时，所研究的各键对相对数目虽然在经过β弛豫后出现了少量的1421，1422和1551键对，但在随后的α弛豫中仅表现出随时间的振荡变化，这说明团簇内原子在这个较高急冷温度下仍表现为如右侧原子结构图所示的无序状态.随着急冷温度降低到600K，虽然团簇内1421和1422键对数目没有发生大的改变，但值得注意的是在进入α弛豫区后，1551键对数目随时间出现增加的趋势，这说明随着急冷温度的降低，一定数量的二十面体局域结构开始出现，这种趋势随着急冷温度的进一步降低表现得更为明显.如图4(c)所示，当急冷温度为200K时，在0.1 ps后，团簇内1551键对数目明显增加，并且这种结构在50 ps后仍有相当数量的存在，这说明大部分的二十面体局域结构可以较长时间地存在.这种1551键对的出现也可以由图右侧100 ps时团簇原子排布出现环状形状得到验证.这里还需指出的是，在急冷温度200K为600K时，1421和1422键对都保持很低的值，这说明团簇内成核的结构主要为二十面体结构，并且这些成核的二十面体没有发生明显的生长.

4. 结论

应用基于嵌入原子方法的正则系综分子动力学研究熔融Cu57团簇冷却到不同温度时随时间变化函数及团簇局域结构的变化.模拟结果表明，尽管熔融Cu57团簇在被急冷到三个不同温度时都出现相似的短时原子集体剧烈运动，但在随后进入α弛豫后，团簇内原子运动及其局域结构随急冷温度不同出现了很大的差异.当急冷温度较高时，α弛豫区的原子运动先是以原子扩散运动为主，然后出现非扩散性的原子局域结构重排.随着急冷温度的降低，在α弛豫前期出现非扩散性的原子局域结构调整，随后原子局域结构的变化来自于原子的扩散性运动，这时出现一定数量的二十面体结构.当急冷温度很低时，在α弛豫区逐渐表现为非扩散性原子局域结构重排，相当数量的成核二十面体结构形成，其中大部分可以保持其局域稳定结构，这些成核的二十面体结构没有发生明显的生长.

［1］Quan H J，Gong X G 2000 Chin.Phys.9 656

［2］Li T X，Ji Y L，Yu S W，Wang G H 2000 Solid State Commun. 116 547

［3］Xu S N，Zhang L，Zhang C B，Qi Y 2007 Acta Metall.Sin.43 379(in Chinese)［徐送宁、张林、张彩碚、祁阳2007金属学报43 379］

［4］Wu H，Desai S R，Wang L S 1996 Phys.Rev.Lett.77 2436

［5］Liu H B，Ascencio J A，Alvarez M P，Yacaman M J 2001 Surf. Sci.491 88

［7］Leutheusser E 1984 Phys.Rev.A 29 2765

［8］Kob W，Andersen H C 1995 Phys.Rev.B 22 4

［9］Kob W，Andersen H C 1995 Phys.Rev.E 51 5

［10］Tokuyama M 2006 Physica A 23 62

［11］Pang H，Jin Z H，Lu K 2003 Phys.Rev.B 67 094113

［12］Kob W，Donati C，Plimpton S J，Poole P H，Glotzer S 1997 Phys.Rev.Lett.79 15

［13］Donati C，Glotzer S C，Poole P H，Kob W，Plimpton S J 1999 Phys.Rev.E 60 3

［14］Gleim T，Kob W，Binder K 1998 Phys.Rev.Lett.81 20

［15］Zhang L，Zhang C B，Qi Y 2007 Chin.Phys.16 77

［16］Yang Q W，Zhu R Z，Wen Y H 2005 Acta Phys.Sin.54 89(in Chinese)［杨全文、朱如曾、文玉华2005物理学报54 89］

［17］Zhang L，Zhang C B，Qi Y 2008 Phys.Lett.A 372 2874

［18］Zhang L，Zhang C B，Qi Y 2009 Physica B 404 205

［19］Zhang L，Sun H X 2009 Solid State Commun.149 1722

［20］Zhang L，Sun H X 2009 Chin.J.Chem.Phys.22 69

［21］Mei J，Davenport J W，Fernado G W 1991 Phys.Rev.B 43 4653

［22］Honeycutt J D，Andersen H C 1987 J.Phys.Chem.91 4950

［23］Clarke A S，Jonsson H 1993 Phys.Rev.E 47 3975

PACC:3640B，7115Q，6185，6800

*Project supported by the State Key Development Program for Basic Research of China(Grant No.G2006CB605103).

†Corresponding author.E-mail:zhanglin@imp.neu.edu.cn

Molecular dynamics study of relaxation and local structure changes in a rapidly quenched molten Cu57cluster*

Fan Qin-Na Li Wei Zhang Lin†
(College of Science，Northeastern University，Shenyang110004，China)
(Received 11 May 2009;revised manuscript received 3 July 2009)

Relaxation and local structure changes of a molten Cu57cluster during rapidly quenching have been studied by molecular dynamics simulation using embedded atom method.With decreasing quenching temperature，atom motion details are analyzed using three factors，including the mean square displacement，incoherent intermediate scattering function，and non-Gaussian parameter，while the local structure changes are identified by pair analysis.Simulation results reveal that after a drastic collective motion of atoms，the temperature greatly affects the relaxation processes of the cooled cluster.At a high quenching temperature，after atoms dramatically move in a β relaxation region，diffusion motion of the atoms plays a dominant roles followed by non-diffusion rearrangements of local atomic structures，and no nucleation occurs.When the temperature decreases，local structure changes of atoms occur as the initial dramatic motion，then through the diffusion of atoms in the α relaxation region，and some unstable icosahedral structures are observed.At a low quenching temperature，the structure changes in the α relaxation region result mainly from non-diffusion rearrangement of the atom positions，and a notable amount of icosahedral structures are formed.

cluster，molecular dynamics，computer simulation，surface

book=102,ebook=102

*国家重点基础研究发展计划(批准号:G2006CB605103)资助的课题.

†通讯联系人.E-mail:zhanglin@imp.neu.edu.cn