merge-fdf.py：高效合并Siesta fdf输入文件的小工具

干这行的人肯定懂：一个siesta计算没跑起来，多半不是物理模型出了问题，而是*.fdf输入文件写得不对。fdf 这东西看着就是个文本，但参数、结构、基组、k 点全混在一个文件里，手工维护起来非常烦。尤其是结构优化完要更新坐标、批量准备掺杂构型、或者想把参数模板和结构模板分开管理的时候，复制粘贴 fdf 段落就成了最磨人的重复劳动。

merge-fdf.py就是为解决这个问题写的一个小工具：按照 fdf 自身的语法规则，把多个 fdf 文件按顺序拼接成一个，同名的普通参数和%block结构块以后出现的为准，位置保持在第一次出现的地方。脚本只用 Python 标准库，不依赖第三方包，能跑python3就能用。适合刚接触 siesta、还在手工拼 fdf 的人，也适合已经在管理大量模板文件、想把输入文件生成流程化的老手。

这篇文章我会把 fdf 的格式拆开讲清楚，然后结合merge-fdf.py的完整实现，说说脚本背后的设计思路，再用三个真实场景演示怎么用。文末还会把我调试过程中踩过的坑和后续扩展方向一并写出来，希望对你有帮助。

1. 为什么需要 merge-fdf.py：三个常见的手拼现场

1.1 结构优化后的坐标替换：参数没变，坐标全变了

跑完一轮结构优化，下一步通常是把优化后的原子坐标拿回来，套用同一套计算参数继续算电子结构、能带或者态密度。这个操作听着简单，实际做起来很折腾。

我早期是打开输出文件，从里面找到最后一步的坐标块，然后手动复制替换到输入 fdf 里。坐标少的时候还好，几十个原子也能凑合，可一旦体系大了，或者做了可变晶胞的优化需要同时更新格矢，手动操作很容易出错。复制的时候漏一行、多复制一个原子、或者把旧坐标的某几列留下来，这类问题几乎没有预兆，直到计算跑起来才发现几何构型已经乱了。

更麻烦的是坐标块替换完之后，还得检查NumberOfAtoms、NumberOfSpecies，以及%block ChemicalSpeciesLabel是否同步。原子数不变还好，假如你在结构文件里多删了一个原子，NumberOfAtoms没同步更新，siesta 通常会直接报错或者读入一个错误的结构。这种低级错误很消磨耐心。

1.2 批量准备掺杂/吸附构型：重复劳动最磨人

做表面吸附、掺杂筛选这类工作时，经常要准备大量结构只有微小差异的输入文件。比如一个表面模型，吸附位点有顶位、桥位、空位三种，每种又分不同吸附分子朝向，算下来几十个目录。每个目录里的 fdf 参数部分完全一样，结构部分只在吸附原子坐标上有区别。

如果靠手工复制，就要先复制一份模板，再改坐标，再另存为下一个目录。重复几十次之后，人很容易麻木，一麻木就会出错。我见过有同事把 A 位点的坐标复制到 B 位点的目录里，跑完一个数据点才发现吸附构型复制错了，整个计算白跑。这不是细心不细心的问题，而是这种流程设计本身就有问题——一旦涉及批量生成输入文件，就必须让脚本代替手工去合并、去覆盖、去生成。

1.3 模板化参数管理：单文件维护最优解

还有一种情况不是批量计算，而是你想把参数管理得更清晰。一个完整的 fdf 里，计算参数、基组定义、k 点设置、结构坐标全都混在一起，想单独调整某个环节就得在几百行里找位置。

更合理的做法是把 fdf 拆成几个逻辑文件，比如params.fdf只管泛函、截断能、收敛判据，species.fdf只管ChemicalSpeciesLabel和 PAO 基组，structure.fdf只管格矢和坐标。需要调整哪部分就改哪个文件，最后用脚本合并成最终输入。这样每个文件都很短，可读性高，也方便 git 做版本管理。

但拆分之后必须有一个可靠的合并工具，否则从"一个文件"变成"五个文件"，反倒给每次提交作业增加了步骤。merge-fdf.py的核心价值就在这里：让文件拆分成为一种日常习惯，而不是负担。

2. 动手之前先拆开 fdf：它到底是个什么格式

2.1 fdf 的两种基本语法：普通键值行和 %block 结构块

要写对拼接逻辑，第一步是搞清楚 fdf 文件里到底有哪几种语法。Siesta 的 fdf 格式大体上就两类，一类是普通键值行，一类是%block结构块。

普通键值行是最常见的形式，一行里包含一个关键词和它的值：

SystemLabel siesta_run MeshCutoff 400 Ry PAO.EnergyShift 25 meV MD.TypeOfRun CG MD.NumCGsteps 200 MD.MaxForceTol 0.01 eV/Ang

这类行遵循"关键词 + 空格 + 值"的规则，值里可能带单位。关键词本身大小写不敏感，写meshcutoff和写MeshCutoff效果一样，但为了可读性通常会保留惯例写法。

%block结构块则用于多行数据，典型的包括：

%block LatticeVectors 3.84000 0.00000 0.00000 0.00000 3.84000 0.00000 0.00000 0.00000 3.84000 %endblock LatticeVectors %block AtomicCoordinatesAndAtomicSpecies 0.00000 0.00000 0.00000 1 1.92000 1.92000 1.92000 2 %endblock AtomicCoordinatesAndAtomicSpecies

再比如%block ChemicalSpeciesLabel、%block kgrid_Monkhorst_Pack（老版本写法）、%block PAO.Basis这些，都是同一种结构块语法。它们的共同特点是：以%block 名字开头，以%endblock结尾，中间是多行内容。这里必须注意，%endblock后面一般可以跟名字，也可以不跟，但写全名字会让文件更清晰，也方便脚本做校验。

fdf 格式的核心语法我用一张表总结一下：

2.2 拼接最容易出错的两个地方

第一是 block 的多行属性。很多人第一次写合并脚本时，会天真地按"行"去拼接：把两个文件的每一行按顺序堆到一起。这在小样本下可能碰巧能跑通，但一旦遇到%block结构块，按行拼接就会把块内容拆得七零八落，甚至把一个文件的后半个 block 和另一个文件的 block 残留拼在一起。真正稳妥的做法是先把文件切成"块"级别的内容单元，再对单元做合并。

第二是结构定义的上下文关系。坐标块%block AtomicCoordinatesAndAtomicSpecies里的原子序号，必须和%block ChemicalSpeciesLabel里定义的物种序号一一对应。如果你合并时只更新了坐标块，却忘了同步物种定义，计算会直接报错。同样，如果你的格矢可变，优化后的LatticeVectors也要一并替换进去，否则会出现坐标和格矢不匹配的物理问题。

2.3 一个关键认知：fdf 的"后者覆盖前者"规则

fdf 文件在解析的时候，对于同一个关键词，多次出现的处理逻辑一般是后面出现的覆盖前面出现的。这一点很关键，因为它直接影响合并脚本的设计方向。

如果你的某个文件里可能出现同一个参数的多个定义，最安全的合并策略就是顺着文件顺序逐个处理，后遇到的同名参数替换之前的同名参数。这样既能保证"无脑拼接"时可能产生的重复定义问题被消除，又能很好地适配"用后面的结构文件覆盖前面的模板文件"这种典型用法。

因此，merge-fdf.py的核心语义可以定成：合并后同一个普通键或同一个 block 只保留一个，以所有输入文件中最后一个出现的为准。

3. merge-fdf.py 的核心实现：先把文件切成"块"，再按覆盖规则组装

3.1 整体设计思路：解析、合并、写出三步走

脚本整体分为三个层次：解析层、合并层、输出层。

解析层负责把单个 fdf 文件解析成一个有序的条目列表。每个条目可能是普通键值行、%block结构块、或者注释/空行。解析层不关心数值是什么，只负责把结构完整地切分出来。

合并层遍历所有输入文件，维护一个"当前合并结果"列表。遇到普通键或 block 时，先看是否已经出现过同名条目，如果出现过就在原位置替换内容，如果没出现过就追加到尾部。注释和空行保持原样追加，保证文件读起来不至于太奇怪。

输出层把合并结果写回文本文件，统一使用 Unix 换行符，避免 Windows 和 Linux 混用时的格式问题。

3.2 fdf 解析器：识别 block 和普通键

解析器的核心是一个顺序扫描循环，用正则判断当前行是%block开始、%endblock结束、还是普通键值行。这里我把普通键值的单位部分也一并保留在原始行里，不拆分，这样后续写出去的时候不会因为重新拼接而丢失单位。

import re import sys def parse_fdf_file(filepath): """把单个 fdf 文件解析成条目列表。 每个条目是一个元组 (kind, name, content, lineno) kind 取值： 'block' -> %block 结构块 'key' -> 普通键值行 'raw' -> 注释或空行 content 保留原始文本内容，方便原样输出。 """ entries = [] with open(filepath, 'r', encoding='utf-8', errors='replace') as fh: lines = fh.readlines() i = 0 n = len(lines) while i < n: line = lines[i] stripped = line.strip() # 空行和注释作为 raw 条目保留 if not stripped or stripped.startswith('#'): entries.append(('raw', '', line, i + 1)) i += 1 continue # 匹配 %block 开头 block_match = re.match(r'^\s*%block\s+(\S+)\s*$', line, re.IGNORECASE) if block_match: block_name = block_match.group(1).lower() buf = [line] i += 1 closed = False while i < n: buf.append(lines[i]) if re.match(r'^\s*%endblock', lines[i], re.IGNORECASE): closed = True i += 1 break i += 1 if not closed: print(f"警告: {filepath} 中的 %block {block_name} 没有闭合", file=sys.stderr) entries.append(('block', block_name, ''.join(buf), 'block')) continue # 普通键值行：第一个空白字符之前的字段当作 key parts = stripped.split(None, 1) if len(parts) >= 1: key = parts[0].lower() entries.append(('key', key, line, i + 1)) else: entries.append(('raw', '', line, i + 1)) i += 1 return entries

注意几个细节。

第一，所有 key 和 block 名都统一转成了小写。因为 fdf 关键词本身大小写不敏感，用小写做指纹可以避免LatticeVectors和latticevectors被当成两个不同条目。

第二，遇到%block时用循环把整个块收集起来，而不是只记一行。这样后续替换时才能保持块内容的完整性。

第三，对于未闭合的 block，脚本会打出一条警告。这个分支看起来不起眼，但实际很有用。我调试的时候遇到过 fdf 文件被截断的情况，如果没有这个发现，合并出来的文件会让 siesta 在解析阶段直接崩溃，排查起来很痛苦。

3.3 覆盖合并逻辑：位置不变，内容替换

解析之后，合并的核心就是一个有"记忆"的顺序遍历。我用一个列表保存最终条目，用一个字典记录每个 key/block 在列表中的位置。这样遍历到同名的新条目时，可以直接替换掉旧条目，同时保持旧条目的原始位置。

为什么保持原始位置很重要？因为 fdf 文件里有些参数之间存在上下文习惯。比如AtomicCoordinatesFormat Ang通常写在坐标块附近，如果你把所有结构块都挪到文件末尾，虽然大部分情况下 siesta 也能读，但读文件的人会觉得很奇怪。保持位置不变，实际效果就是"原位更新"，和手动编辑文件的直觉一致。

def merge_files(input_files): """按顺序合并多个 fdf 文件，返回合并后的条目列表。""" merged = [] index = {} for filepath in input_files: for ent in parse_fdf_file(filepath): kind, name, content, _ = ent[0], ent[1], ent[2], ent[3] if kind == 'raw': # 注释、空行直接保留 merged.append(ent) continue # block 和 key 用不同的前缀做指纹，防止互相覆盖 fingerprint = ('b:' if kind == 'block' else 'k:') + name if fingerprint in index: pos = index[fingerprint] merged[pos] = ent # 替换旧条目，位置不变 else: index[fingerprint] = len(merged) merged.append(ent) return merged def write_fdf(entries, out_path): """把条目列表原样写回文件，统一使用 Unix 换行。""" with open(out_path, 'w', encoding='utf-8', newline='\n') as fh: for ent in entries: fh.write(ent[2])

这十几行就是合并逻辑的全部。它没有做任何数值层面上的处理，因为拼接 fdf 的底线是不改变原文的含义，只做结构上的重组。

3.4 命令行入口：一个文件搞定，不依赖第三方库

为了让脚本在命令行下用起来顺手，我再加一个简单的参数入口。支持指定多个输入文件和一个输出文件，并增加一个--check选项来做基础校验。

import argparse def check_consistency(entries): """检查合并结果中 NumberOfAtoms 与坐标块行数是否一致。""" natoms = None coord_lines = None for kind, name, content, _ in entries: if kind == 'key' and name == 'numberofatoms': parts = content.split(None, 1) try: natoms = int(parts[1]) except (IndexError, ValueError): natoms = None if kind == 'block' and name == 'atomiccoordinatesandatomicspecies': lines = [ln for ln in content.splitlines() if ln.strip() and not ln.strip().startswith('%')] coord_lines = len(lines) if natoms is not None and coord_lines is not None and natoms != coord_lines: print(f"警告: NumberOfAtoms 为 {natoms}，坐标块却有 {coord_lines} 行", file=sys.stderr) return False return True def main(): parser = argparse.ArgumentParser( description='merge-fdf.py: 合并多个 Siesta fdf 文件' ) parser.add_argument('inputs', nargs='+', help='输入 fdf 文件，按从左到右的顺序合并') parser.add_argument('-o', '--output', default='merged.fdf', help='输出文件路径，默认 merged.fdf') parser.add_argument('--check', action='store_true', help='合并后检查原子数与坐标行数是否一致') args = parser.parse_args() entries = merge_files(args.inputs) if args.check: check_consistency(entries) write_fdf(entries, args.output) print(f"已生成 {args.output}") if __name__ == '__main__': main()

至此，merge-fdf.py从解析到输出就完整了。整份代码不到一百行，没有什么高深算法，但它在实际工作流里的用处非常大。

4. 实战演示：三种场景下 merge-fdf.py 的具体用法

4.1 场景一：结构优化后更新坐标

假设你有一个base.fdf，里面是完整的计算参数和初始结构。优化结束后，你想保留所有参数，只把坐标和格矢换成优化后的最终结构。

这时你可以先把最后一帧结构单独抽成一个final_struct.fdf。比较省事的方式是用 ase 读取 Siesta 的轨迹文件，把最后一帧原子坐标写出来，再拼成一个干净的结构块。下面是示例代码：

from ase.io import read atoms = read('siesta.MD', index=-1) with open('final_struct.fdf', 'w') as f: f.write('AtomicCoordinatesFormat Ang\n') f.write('%block AtomicCoordinatesAndAtomicSpecies\n') for atom in atoms: # 注意: species 编号必须与 ChemicalSpeciesLabel 一致 sp = atom.tag if atom.tag else 1 pos = atom.position f.write('%20.10f %20.10f %20.10f %d\n' % (pos[0], pos[1], pos[2], sp)) f.write('%endblock AtomicCoordinatesAndAtomicSpecies\n')

如果你的优化过程允许晶胞变化，还需要把最后的LatticeVectors也提取出来。很多版本的siesta.MD里会周期性打印格矢，提取时注意取最后一帧对应的格矢，不要取初始格矢。

然后执行：

python3 merge-fdf.py -o run.fdf base.fdf final_struct.fdf

这样生成的run.fdf，参数部分完全继承自base.fdf，坐标部分自动替换成最后一帧结构，所有结构块的原始位置保持不变。执行完之后建议加一个--check：

python3 merge-fdf.py -o run.fdf base.fdf final_struct.fdf --check

这个检查会统计NumberOfAtoms和坐标块行数，万一你在提取结构时漏了一个原子，它会在提交作业之前就提醒你。

4.2 场景二：批量生成吸附构型输入文件

批量场景下，merge 的优势更加明显。你只需要维护好两个文件夹：一个是存放公共参数模板的目录，一个是存放各构型结构文件的目录。

比如结构文件命名为site1.fdf、site2.fdf、site3.fdf，每个文件里只写坐标块和格矢块：

for site in site1 site2 site3; do python3 merge-fdf.py -o run_${site}.fdf params.fdf ${site}.fdf done

循环跑完，run_site1.fdf、run_site2.fdf、run_site3.fdf就都生成好了。每个文件里参数部分完全一样，结构部分来自对应的构型文件。接下来不管是本地批量跑还是丢到集群，都可以接一个作业生成脚本继续往下走。

这里有一个需要注意的地方：如果某个构型文件里没有写NumberOfAtoms，那么合并后的文件会沿用params.fdf里的NumberOfAtoms。如果你的几个构型原子数不同，就必须在各自的结构文件里显式写入NumberOfAtoms，否则合并结果就是错的。这个容易出现，我建议在每个结构文件头部都写上原子数和物种数，不要依赖模板去猜。

4.3 场景三：参数模板与结构模板分离管理

我现在的日常工作方式，是给每个项目建一个templates/目录，里面放几类模板：

需要生成一个新作业时，按需合并：

python3 merge-fdf.py -o newjob.fdf params.fdf species.fdf kpoints.fdf cell.fdf

这样做的收益是，参数调整只动params.fdf一个文件，结构替换只动cell.fdf一个文件，其他人接手时也能一眼看清整个计算设置的分层结构。配合 git 之后，每次计算的输入文件从哪来、改过什么，都清清楚楚。

5. 实测踩坑：从"能跑"到"跑对"

5.1 编码和换行符：Windows 上打开 Linux 文件容易翻车

siesta 主要跑在 Linux 集群上，但很多时候模板文件是在 Windows 上编辑的。Windows 的文本编辑器和 Linux 的换行符不一致，文件里会混入\r\n。如果不做处理，merge 生成的 fdf 可能在某些环境下解析出奇怪问题。

我在write_fdf里显式指定了newline='\n'，读取时用了encoding='utf-8', errors='replace'。读取时用errors='replace'可以避免遇到非法编码时直接崩溃。这个小细节带来的好处是，无论输入文件是 Windows 还是 Linux 格式，输出文件统一是 Unix 换行，提交到集群上不会因为换行符问题报错。

需要说明的是，我并没有在读取时把\r\n显式转成\n，因为 Python 打开文本文件时会做通用的换行转换，只要write_fdf输出时统一成\n就够了。实际用下来这个方案最省心。

5.2 坐标精度不能顺手丢

写解析器时，一个很容易犯的错误是对坐标数值做格式化。比如读取坐标块里的每一行，把字符串转成 float 再拼回去。浮点数在打印时如果位数不够，会悄悄丢失精度。对原子位置来说，哪怕丢掉一位有效数字，优化可能就白做了。

所以我一开始就没打算在解析阶段去理解数值的内容。整个脚本把普通键值和 block 的内容都当作"透明文本"处理，原样保留、原样输出。坐标位数在用户自己的结构文件里是多少，合并出来就是多少，脚本绝不碰数值。这是一个很重要的设计取舍：拼接脚本只负责结构和去重，不负责格式化数据。

5.3 同名键在不同位置覆盖的陷阱

后出现的同名条目会替换先出现的同名条目，这个规则大多数时候很好用，但也会带来一个隐蔽的问题。比如params.fdf里你已经写好了SystemLabel run_001，后面的kpoints.fdf里不小心也写了一个SystemLabel run_kpt，合并之后输出文件里的系统标签就变成了run_kpt，作业输出前缀全变了。

这种问题不算致命，但会浪费一轮排查时间。我的对策是两层：第一，所有模板文件保持简洁，最上层的公共参数只在一个文件里出现；第二，在用脚本生成最终输入后，打开文件扫一眼开头和结尾，确认没有明显不符合预期的重复定义。如果文件很多，也可以用grep -i systemlabel快速确认。

5.4 合并后必做的两道检查

第一道检查是差异对比。可以把合并结果和原始的 base 文件做一次 diff，看看到底哪些行被替换了、哪些行被追加了。特别是在更新结构坐标时，diff 能很明显地暴露坐标块是否被正确替换，以及NumberOfAtoms是否出现了重复定义。

第二道检查就是前面提到的--check。这个检查不复杂，但对吸附构型这类原子数可能变化的场景非常管用。我建议在任何更新过结构的合并操作中都加上--check，成本几乎为零，却能提前发现最典型的低级错误。

6. 扩展思路：让它长成你自己的小工具

6.1 给脚本加一个 --diff 模式

实际用起来之后，你会发现"合并前看差异"和"合并后看差异"是两种完全不同的需求。有些模板文件很长，你只想知道两个版本之间到底改了什么。

可以在脚本里增加一个--diff参数，直接调用系统diff命令对比两个合并结果。更简单的方式其实是把merge-fdf.py和命令行 diff 组合使用：

python3 merge-fdf.py -o old.fdf params_old.fdf cell.fdf python3 merge-fdf.py -o new.fdf params_new.fdf cell.fdf diff old.fdf new.fdf

这样能快速看到参数变更对整体输入的影响，适合在做收敛性测试时使用。

6.2 反向需求：从一个 fdf 里单独抽出结构块

拼接有了，反向提取也是常用功能。有时候你从别人的项目里拿到一个完整的 fdf，想把里面的结构部分抽出来作为自己的模板，配合参数文件重新组织项目。

这个需求不需要再写一个独立工具，只要用parse_fdf_file把 block 挑出来：

from merge_fdf import parse_fdf_file for kind, name, content, _ in parse_fdf_file('example.fdf'): if kind == 'block' and name in ('latticevectors', 'atomiccoordinatesandatomicspecies'): print(content)

有了解析层，这类小需求都可以顺手实现，不用再从零写正则。

6.3 与自动化作业生成流水线结合

最后说一个我个人的工作流习惯。我现在不会手动运行merge-fdf.py，而是把它嵌到整个作业生成流水线里。大体流程是：先用 Python 脚本批量生成所有结构文件，再用循环调用merge-fdf.py把每个结构的 fdf 拼好，接着生成对应的 PBS/SLURM 作业脚本，最后统一提交。

这样的好处是，每个目录下的输入文件都可以随时从模板重新生成。即使某个目录被误删了，只要跑一遍生成脚本，所有输入就能恢复原样。对我来说，merge-fdf.py不再是一个孤立的脚本，而是整个计算流程中一个可靠的基础设施。

如果你也经常和 fdf 文件打交道，建议从最基础的合并开始用，慢慢加上适合自己项目的参数和扩展。这个脚本不大，但它能把那些重复、容易出错的手工操作从日常工作中彻底去掉，省下来的时间足够你多跑好几个构型了。