CO在Pt(111)表面吸附的VASP计算全流程：从Materials Studio建模到吸附能分析-平芜编程栈

CO在Pt(111)表面的吸附，几乎是每个刚开始做多相催化计算的人都要过一遍的经典体系。体系本身看着简单——一个表面、一个CO分子、一次结构优化、一个吸附能——但真正动手做起来，从Materials Studio建模到VASP参数设置，中间的细节和坑一点都不少。这个系列的第一篇，我就把完整流程拆开讲清楚：为什么用这个模型、怎么在MS里搭出Pt(111)表面、怎么摆CO分子、VASP输入文件怎么配、吸附能怎么算，以及我踩过的那些坑。

这篇文章适合谁？刚接触VASP的研究生、想自学表面吸附计算的本科生，以及想快速复现CO/Pt体系结果的科研人员。如果你已经能熟练跑通简单的VASP任务，也可以跳过基础部分，直接看第5章的能量解读和第6章的排错清单。下面我们直接开始。

1. 为什么是CO/Pt(111)：从催化问题到计算模型

1.1 这个体系在催化里到底有多重要

CO在铂表面的吸附，是一切铂基催化剂故事的开端。汽车尾气三元催化器里，CO氧化是核心反应；质子交换膜燃料电池的铂电极上，微量CO就能让催化剂“中毒”失活；水煤气变换、甲醇氧化这些工业反应，第一步几乎都避不开CO的吸附活化。可以说，理解了CO怎么吸附在Pt表面，就相当于拿到了铂基催化表面化学的入门钥匙。

在计算化学领域，CO/Pt(111)的地位则更特殊：它是表面科学里被研究得最充分的标准模板体系。从上世纪单晶表面的实验测量，到后来第一性原理计算的广泛benchmark，无数篇方法学论文拿它做测试体系。对新入坑的计算研究者来说，这个体系足够“干净”：吸附分子小、表面结构稳定、吸附位点有限，非常适合用来建立对吸附计算全流程的直觉。我当年入门时，导师让我做的第一个任务就是复现CO在Pt(111)四个位点的吸附能，今天回头看，这个选择确实很值。

1.2 为什么选Pt(111)而不是别的表面

铂是面心立方（FCC）结构，常见低指数面有(100)、(110)、(111)。实验上，Pt(111)是热力学最稳定的晶面，也是最容易制备出高质量单晶的表面，所以大量表面科学研究都集中在它上面。理论计算选它还有一个实际原因：FCC(111)面原子排列紧凑、对称性高，一个比较小的超胞就能描述表面周期性，计算量控制得住。

这里顺便解释一个容易混淆的概念：Pt(111)这个晶面指数指的是表面法线方向沿着晶体的[111]方向，FCC晶体里这个方向是原子密排面。密排面上的原子配位数为9（体相是12），表面悬键少，所以它比(110)和(100)更稳定。这也意味着，理论计算里用(111)面得出来的吸附趋势，通常和真实多晶催化剂上占主导的晶面行为最接近。

1.3 整体计算流程设计

整个项目可以拆成六步闭环，每一步都依赖前一步的输出，别跳步：

用Materials Studio构建Pt晶胞，切出(111)表面，加真空层，构建超胞。
在表面上放置CO分子，生成top、bridge、fcc、hcp四种初始吸附构型。
准备VASP四个输入文件：POSCAR、POTCAR、INCAR、KPOINTS。
运行结构优化，让原子位置充分弛豫。
分别计算干净表面和孤立CO分子的参考能量。
按公式计算吸附能、比较位点稳定性、检查结构合理性。

这套流程看似简单，但每一步都有细节。下面我按顺序拆解，重点写建模时“为什么这么设”的逻辑。

2. Materials Studio建模：从晶胞到Pt(111)表面

2.1 获取Pt晶体结构和切面操作

Materials Studio（MS）里建Pt表面，最省事的方式是直接从结构库导入体相结构。具体路径是File -> Import，在Structures/metals目录下找到Pt，或者用菜单里的Build -> Crystals -> Build Crystal手动输入晶格常数。Pt的FCC晶格常数实验值约3.92 Å，这个值可以直接作为初始结构使用；如果你后面做高精度计算，更建议先用VASP优化一遍体相结构，再拿优化后的晶格常数去切表面，这样能消除晶格常数的误差影响。

导入体相Pt后，关键操作是切面。菜单路径：Build -> Surfaces -> Cleave Surface。在弹出的面板里，把Cleave Plane的Miller指数改成(1, 1, 1)，Surface mesh和Top thickness按需设置。Top thickness的物理含义是slab厚度，我建议先设4层（约7.1 Å）。点击Cleave之后，屏幕上会出现一个重复的表面结构，注意观察坐标系：MS会自动把a轴和b轴放在表面平面内，c轴垂直于表面。如果你的MS界面左下角有晶体轴指示器，务必确认c轴方向是垂直向上的，这是后面加真空层和导出POSCAR不出错的前提。

切面完成后，Build -> Symmetry -> Supercell可以建超胞。我建议的默认方案是2x2超胞，后面会详细解释为什么。

2.2 真空层与超胞尺寸的设置逻辑

切出来的表面是一个周期性延拓的slab，沿c方向会无限重复。如果不加真空层，相邻slab之间会直接“粘”在一起，这显然不符合真实表面条件。我们需要在slab上方留出一段真空区域来隔离周期性镜像。操作是Build -> Crystals -> Build Vacuum Slab，真空层厚度一般设15 Å。

为什么是15 Å而不是5 Å或者30 Å？真空层太薄，slab上表面的电子密度会和周期性相邻的slab下表面发生重叠，带来虚假的相互作用能；太厚则纯粹浪费计算资源，因为平面波基组会把真空区域也一起描述，需要的平面波数量会增加。15 Å是一个非常成熟的折中值，在大多数赝势下都能让静电相互作用和电子密度在中间区域衰减到可忽略程度。如果你计算带偶极的吸附体系，可能需要更厚的真空层或者配合偶极修正，这个后面第5章会提。

超胞尺寸决定的是吸附覆盖度。2x2超胞意味着原胞里有4个表面Pt原子，吸附1个CO分子时覆盖度就是1/4 ML；3x3超胞则是1/9 ML。覆盖度越高，相邻CO分子之间的排斥作用越强，吸附能绝对值会变小。实验上研究单分子吸附行为，倾向于低覆盖度，但计算时要平衡精度和成本。我建议入门阶段用(2x2)，理由有三：一是这个尺寸下CO—CO最近间距约5.5 Å，排斥作用已经不算太强；二是4个表面原子的模型跑起来很快，新手迭代方便；三是绝大多数文献都有(2x2)数据可以直接对比。如果后续你要发表文章或对比实验低覆盖度数据，可以再升级到(3x3)。

2.3 表面原子层数与固定策略

前面说切4层，这里展开讲层数的影响。slab层数太少，上下两个表面之间会通过金属内部电子态发生耦合，这种“量子尺寸效应”会导致表面能和吸附能出现非物理的波动；层数太多又线性增加计算成本。对于Pt(111)，4-5层已经能给出收敛的结果。如果你担心层数影响，可以做一个简单的收敛性测试：比较4层、5层、6层的CO吸附能，看差值是否小于20 meV。实测下来，Pt(111)在4层到5层之间通常就能收敛得很好。

固定策略是很多新手容易忽略的点。我们在MS里切出的slab其实是一个上下都有表面的孤立薄膜。真实晶体中，表面下层还有无限延伸的体相结构，它们会对表面施加强约束，防止整个slab在优化时发生整体平移、弯曲或滑动。为了模拟这个约束，惯例是固定slab底部一定厚度的原子层，只放开顶层和吸附分子弛豫。对于4层slab，通常固定最下面1-2层；对于5-6层slab，固定一半左右也算常见。

在MS里固定原子，可以直接选中底部原子（框选或按原子序号选），然后在属性面板里找到约束选项，勾选固定XYZ坐标；也可以在导出POSCAR后，用文本编辑器修改Selective dynamics标志，把要固定的原子行改成 F F F，其余原子 T T T。我更推荐后者，因为最终VASP认的是POSCAR里的标志，手动检查一遍可以避免MS导出被静默丢掉的意外。

3. 构建CO分子与四种典型吸附位点

3.1 CO分子的构建与取向选择

CO分子很小，在MS里可以用Sketch工具手动画出来，也可以从结构库里找现成分子。画的时候注意两点：一是C—O键长约1.13 Å，不用画得特别精确，结构优化会自己调整；二是分子要放在一个独立的周期性盒子里，这个盒子尺寸最好和后面的吸附超胞保持一致，原因后面第5章讲参考能量时会说。

CO吸附在过渡金属表面，几乎总是以碳端朝下的方式成键。原因要从分子轨道角度理解：CO的HOMO是5σ轨道（孤对电子在C端），LUMO是2π反键轨道；金属表面给电子到CO的2π轨道，同时CO的5σ轨道向金属表面提供电子，这种“σ给电子-π反馈”键合模式，要求C原子指向金属。反过来让O端朝下的“倒吸附”能量会高很多，不属于我们关心的稳定结构。

另外强调一下分子轴取向。初始构型里CO分子轴最好垂直于表面，也就是C在吸附位点正上方、O垂直朝外。为什么？因为垂直构型对称性高、收敛稳定，即使真实最低能量结构有一点倾斜，优化程序也能从垂直起点自然弛豫过去。一上来就故意摆一个角度，反而可能让优化过程多走弯路。

3.2 识别fcc、hcp、top、bridge位点

Pt(111)表面上存在四类高对称吸附位点，这是整个体系的核心知识点：

top（顶位）：CO吸附在单个Pt原子正上方。
bridge（桥位）：CO桥接两个相邻Pt原子，C位于两个Pt连线中点。
fcc hollow（FCC空位）：由三个表面Pt原子围成的三角形空洞，这个空洞正下方没有原子。
hcp hollow（HCP空位）：同样由三个表面Pt原子围成的空洞，但这个空洞正下方正好有一个第二层Pt原子。

fcc和hcp是新手最常搞混的一对。记住一句话：看空洞正下方有没有原子。Pt(111)的FCC堆积中，表面第一层是A层，第二层是B层，第三层是C层。由A层三个相邻原子围成的“向上三角形”空洞，下面正好对着B层原子，这就是hcp位；而“向下三角形”空洞下面没有原子，直接露出C层空位，这就是fcc位。在MS里把显示模式切到球棍模型，旋转到侧视角，一眼就能分辨。

3.3 初始构型的摆放技巧与注意事项

在MS里创建四个构型，本质上是“复制结构-移动分子-保存”的重复劳动。具体操作可以是：先把CO分子移动到接近目标位点的高度，再用坐标微调面板把C原子的分数量（x,y）放到精确位置。不同位点的定位方法：

top位：把C的x,y坐标设成和某个表面Pt完全相同，高度约2.0 Å。
bridge位：选择两个相邻Pt原子，取它们的x,y坐标平均值作为C的位置，高度约1.8 Å。
fcc/hcp位：选择对应三角形空洞的三个Pt原子，取坐标平均值作为C的位置，高度约1.7-1.8 Å。

O原子的初始位置，建议放在C正上方1.15 Å处，这样CO分子轴正好垂直表面。摆放完后，一定要旋转检查一遍：确认C没有和表面Pt重叠，O没有被塞进slab里。我曾经犯过一个错误，把CO整个放反了，O端朝下还跑了一整轮优化，结果吸附能算出来是正值，才发现构型摆错了。

四个构型建议分别保存为top.xsd、bridge.xsd、fcc.xsd、hcp.xsd，千万别放在一个文件里，后面VASP计算目录管理会非常混乱。

4. VASP输入文件与参数设置详解

4.1 四个输入文件：POSCAR、POTCAR、INCAR、KPOINTS

VASP计算一个体系，需要四个输入文件，缺一不可。POSCAR描述晶格和原子坐标，MS里可以直接导出：File -> Export，格式选择POSCAR。导出后打开看一眼，确认里面坐标格式是分数坐标还是笛卡尔坐标，顺序是否与你的预期一致。

POTCAR是赝势文件，它的元素顺序必须与POSCAR中的原子顺序严格对应。这是新手最容易踩的坑：MS导出的POSCAR原子顺序可能是按元素种类自动排的，而POTCAR是你手动cat的，一旦顺序不一致，VASP会按POTCAR的顺序去读取POSCAR的原子，结果就是电子数完全错乱、能量毫无意义。稳妥做法是：先确定POSCAR里每类元素的书写顺序，再按这个顺序依次cat对应元素的POTCAR。例如：

cat pot_POTCAR/Pt/POTCAR pot_POTCAR/C/POTCAR pot_POTCAR/O/POTCAR > POTCAR

INCAR是计算控制参数文件，核心内容我下一小节单独讲。KPOINTS定义布里渊区采样网格，在表面体系里它的设置有一点特殊性：因为c方向是真空层，没有周期性电子色散，所以k点网格在这个方向只取1。一个典型的Gamma-centered KPOINTS文件长这样：

Automatic mesh 0 Gamma 3 3 1 0 0 0

这里3 3 1分别对应a、b、c三个方向的k点数目，c方向为1就是刻意不采样真空方向。k点数目的选择逻辑在第4.2节展开。

4.2 关键参数的选择逻辑：ENCUT、KPOINTS、ISMEAR

INCAR里最关键的一批参数，我逐个说明背后的考虑。

ENCUT：平面波截断能，直接决定基组完备度。VASP推荐的POTCAR里会写出每个元素推荐的截断能ENMAX，你设置的ENCUT不能低于所有元素里最高的ENMAX。一个实用做法是直接取400 eV，对Pt、C、O组成的体系，这个值已经足够收敛。如果要做严谨测试，可以跑450 eV和500 eV的对比，看总能量差是否小于1 meV/原子。注意，吸附能是两个大能量相减，基组误差需要格外小心，前后计算用相同ENCUT是基本底线。

KPOINTS：对于(2x2)的Pt(111)表面，我推荐用3x3x1的Gamma-centered网格起步。理由很简单：金属体系在费米面附近的态密度变化剧烈，需要一定的k点密度才能准确描述；但表面原胞在a、b方向只有两倍关系，3x3对应布里渊区区段已经足够粗筛。你可以在优化完成后做一个5x5x1的静态计算对比能量，如果吸附能差小于10 meV，说明k点收敛。注意要用Gamma-centered网格，因为表面原胞的布里渊区中心有高对称点。

ISMEAR：这个参数控制能带填充的smearing方法，选错了会报错或者能量不准。金属体系推荐ISMEAR=1（Methfessel-Paxton一阶），配合SIGMA=0.2；如果体系含有很大的真空区和孤立分子，可以用ISMEAR=0（Gaussian smearing）。千万不要对金属用ISMEAR=-5（四面体方法），因为四面体方法在金属体系的k点积分中可能不稳定，尤其优化过程中容易出现负占据数警告。

其它重要参数建议：EDIFF=1E-5控制电子步收敛精度，EDIFFG=-0.02控制离子步收敛力标准；IBRION=2（共轭梯度）和ISIF=2（只优化原子位置，保持晶格不变）是结构优化的标准组合；NELM=200适当提高电子步上限，避免金属体系偶然不收敛。如果担心Pt表面出现磁性，可以加上ISPIN=2试算一次，但CO和Pt(111)通常都是非磁基态，实测算下来磁矩基本为零。

4.3 结构优化收敛标准与流程

我建议把结构优化分成两轮，而不是一步到位。

第一轮，优化干净的slab。把从MS导出的不含CO的POSCAR提交一轮结构优化，让表面原子弛豫到平衡位置。这一步的意义有两个：一是检验slab模型和参数设置是否正常；二是得到一个干净的驰豫结构和总能量，后面算吸附能直接用这个能量作参考。优化结束后，把CONTCAR备份，后续如果需要“干净表面参考体系”，就拿这个CONTCAR当POSCAR用。

第二轮，在驰豫好的slab上放置CO分子（所以前面MS建模型时，最好先用干净slab导出POSCAR，再往MS里加CO构型，分别导出），然后把整个吸附体系提交优化。注意INCAR里底层原子依旧需要固定，你可以在第二轮POSCAR里保留Selective dynamics标志，让底层为F F F，表层和CO为T T T。这样优化的物理图景是：体相衬底不动，CO和表面顶层充分驰豫。

怎么看优化收敛了？打开OSZICAR文件，观察离子步循环里dE和EDIFFG对应的力项是否达到预设阈值。VASP还有一个判断技巧：看连续两步的能量变化是否已经非常小。如果发现离子步在震荡（能量在某个值附近来回跳就是不降），先不要慌，把EDIFFG放宽到-0.05，或者改用IBRION=1（RMM-DIIS），往往能打破震荡。

5. 吸附能计算与结果解读

5.1 吸附能公式与参考能量

吸附能是衡量吸附强度的核心指标，定义如下：

E_ad = E(CO/slab) - E(slab) - E(CO)

其中E(CO/slab)是优化后的吸附体系总能量，E(slab)是同一个超胞尺寸下干净slab优化后的总能量，E(CO)是孤立CO分子在相同盒子中的总能量。吸附能为负表示放热吸附，绝对值越大表示结合越强。

这里有两个容易犯错的关键点。

第一，E(CO)必须用和吸附体系相同尺寸的周期性盒子计算。为什么？CO虽然近似孤立分子，但在周期性边界条件下仍然会与自身的周期镜像发生微弱相互作用。如果CO参考计算用的是10 Å的小盒子，而吸附体系用的是15 Å的超胞，盒子尺寸差带来的能量误差就会进入吸附能。最稳妥的做法是：把CO分子放在一个和slab超胞完全一样的格子里（比如(2x2)Pt(111)的晶格常数的abc），只放一个CO，做一次Gamma-only计算。

第二，E(slab)必须是相同k点设置、相同slab层数、相同固定方式下的干净表面能量。固定原子的处理如果不一致，吸附能里就会混入一层“假能量”。所以我在第4.3节强调，干净slab和吸附体系要用同一套参数跑，最好用同一个父目录下的输入文件改出来的。

举个例子方便理解。假设跑完得到三组能量：E(CO/slab) = -153.42 eV，E(slab) = -152.26 eV，E(CO) = -14.18 eV，那么：

E_ad = -153.42 - (-152.26) - (-14.18) = -1.02 eV

这个数值可以直接和文献同条件计算结果比较。注意啊，eV和kJ/mol换算倍率是96.485，所以-1.02 eV对应约-98 kJ/mol。

5.2 四种位点的能量比较与零点能说明

在PBE泛函水平，常见的计算结果排序是：fcc位点最稳定，hcp位点紧随其后（通常高0.02-0.08 eV），bridge位点再高一些，top位点一般最高。吸附能数值通常在-1.6到-1.9 eV之间（(2x2)超胞PBE结果），不同文献因参数设置会略有浮动。

但这里有一个著名的争议点，值得所有做CO/Pt的人了解：PBE泛函会高估CO在Pt(111)上的吸附强度，同时过度偏好fcc位点。实验上在低覆盖度下观察到top位点具有相当重要的布居，而PBE却强烈预言fcc最低。这不是你计算错了，而是GGA泛函对CO的5σ/2π*轨道相互作用的描述存在系统误差。如果文章里需要更准的位点相对能量，可以考虑RPBE泛函、混合泛函，或者加入DFT-D3色散校正，也可以参考实验吸附热做对比。初学者不必过度焦虑这个偏差，理解现象本身已经超过很多人了。

另外要说清楚：吸附能是0 K电子总能差，没有包含零点能（ZPE）修正和振动熵。如果要严格对比实验吸附焓，需要做频率计算拿到振动模式，再计算ZPE修正。一般来说，CO在Pt上的ZPE修正大约在0.05-0.1 eV量级，对位点相对能量的影响虽然不大，但发表高水平文章时应该补充。

5.3 注意事项：BSSE、偶极修正、色散

平面波+PAW方法相比局域基组有一个天然优点：基组几乎不会出现基组重叠误差（BSSE）的问题，所以常规VASP吸附计算不需要做counterpoise修正。但偶尔遇到体系里有较重元素或使用超软赝势时，还是会有人争论BSSE影响，我的建议是：先用默认参数跑，如果审稿人要求可以再补测，不需要在一开始就过度处理。

偶极修正：slab有真空层时，如果吸附分子（如CO）垂直于表面，体系在c方向会形成净偶极，这样周期性镜像之间会有寄生静电作用。为了消除它，可以在INCAR里加IDIPOL=3和DIPOL参数。对CO/Pt这种净偶极不算太大的体系，开了更稳妥，通常对吸附能的影响在几十meV以内。

色散修正：PBE本身不描述长程范德华色散。对CO这类化学吸附为主的体系，色散贡献较小，可以不加；但如果你要研究含长烃链分子、分子晶体或弱吸附体系，建议用DFT-D3（INCAR里设置IVDW=10、11或12，视VASP版本而定）。CO/Pt体系里，加不加D3对fcc位点的吸附能可能有几十meV的影响，但对位点排序趋势影响有限。我入门时直接用PBE，等后面做系统研究再统一加D3。

6. 常见问题与排查实录

6.1 Materials Studio常见故障：主机名解析与许可证连接

Materials Studio在课题组里是典型的Client-Server架构，想正常打开绘图窗口必须先连上Gateway和许可服务器。新手最常见的两个报错，一个是“系统无法解析主机名称”，另一个是“连不上许可服务器”。我见过的情形多半是这么几类：

第一，许可证服务器的主机名拼写错误或大小写不匹配。License设置里填的服务器名必须和实际机器名完全一致，有时候教师机的机器名带着奇怪的后缀，容易看错。排查命令就是先ping服务器主机名，看能不能解析出正确IP。

第二，客户端的hosts文件缺少映射条目。部分组内机器通过局域网固定IP使用License，但如果服务器IP变了而hosts文件没更新，就会出现“无法解析主机名称”。这种情况检查服务器IP、更新hosts文件映射，通常马上能好。

第三，Gateway服务没启动或端口被防火墙拦截。MS的Gateway在启动时会占用特定端口（默认和License端口有关联），改组网络策略后常常会出现端口不通。这时候检查局域网防火墙设置、确认Gateway服务状态，重启Gateway多半能解决。还有一类很常见的情况是电脑休眠导致Gateway进程假死，重启电脑是最高效的办法。

这个系列后续如果需要，可以专门写一篇MS许可证排查的口袋指南。这里先点到为止。

6.2 VASP计算中的典型报错

第一个高频报错：ZBRENT: fatal error in bracketing。看着吓人，其实本质是结构太差，优化找不到能量下界。原因通常是初始构型里的原子重叠或者CO离表面太近。解决办法很简单：把初始高度恢复到1.8-2.0 Å附近，重新优化。

第二个高频问题：电子步不收敛，OSZICAR里的电子步卡在NELM上限。金属体系偶尔会遇到这种情况，建议把NELM提高到200，减小SIGMA到0.1尝试，或者把ALGO从Fast换成Normal。注意SIGMA减小后smearing变窄，电子步会更难收敛，要同时权衡。

第三个高频问题：离子步震荡停不下来。打开OUTCAR或out文件观察，如果原子坐标在某个值附近来回跳，试试把IBRION从2改成1，或者把EDIFFG放宽到-0.05先跑一个“粗优化”，再拿粗优化结构做严格收敛。

第四个算常见但容易被忽略的问题：优化结束后CONTCAR没有备份。我用ASE提交任务时养成一个习惯，每个任务跑完后第一时间把CONTCAR复制成POSCAR，以便后续做静态计算或继续优化。这个习惯帮我避免过不少“重新跑一遍”的悲剧。

6.3 独家避坑清单

POTCAR顺序必须与POSCAR严格一致，这是检查表里排第一的硬指标。
每个位点单独建目录，目录名用top、bridge、fcc、hcp这种无歧义命名，别用final_final_v2这类的名字。
计算前先确认POSCAR里原子数量：干净slab有多少原子，加CO后增加了几个，心里有数再跑。
结构优化结束后，用VESTA或ASE打开CONTCAR看一眼，确认CO没有“飞走”、表面没有乱重构。出现离谱构型时多数是初始结构问题，不是算法问题。
做吸附能参考能量时，slab参考和CO参考都要用相同ENCUT和相同盒子尺寸，这是能用能量差说话的前提。

再说一个我自己的经验：新手第一次跑CO/Pt(111)时，别贪多，先把一个位点打通跑顺，再看其他三个位点。等到四个位点都出结果，你得到的不仅是一份吸附能清单，更是对整个“建模—计算—分析”闭环的肌肉记忆。后面再接触过渡态搜索、反应路径、动态模拟，都是在今天这条路上一层层往上搭。把这一步走扎实，VASP吸附计算的骨架就算立住了。