企业官网建设流程全解析

简介：这份资源围绕第一性原理计算软件Octopus展开，聚焦时间依赖密度泛函理论（TDDFT）模拟计算，面向具备一定量子化学与凝聚态物理基础、希望上手电子结构模拟的科研人员与研究生。内容涵盖基态DFT预计算、TDDFT参数设置、激发态响应计算、吸收与荧光光谱分析以及电子密度与轨道可视化等完整流程，可帮助读者理解分子与材料光学性质的模拟思路。压缩包为octopus-4.1.0版本，整体约4.58MB，通常包含源代码、编译脚本、文档与示例输入文件等，便于自行编译配置计算环境。目前已有409人学习下载，适合作为TDDFT入门与实操的参考材料，配合官方文档逐步掌握参数设置与输出解读，为太阳能电池材料、光催化剂等方向的研究提供计算支持。

1. octopus 与 TDDFT 模拟计算：从一束激光到电子云的真实响应

很多做第一性原理的同行第一次接触含时密度泛函理论（TDDFT）时，都会卡在同一个地方：基态算得挺顺，一到“加光场、看激发、追电子动力学”就翻车。要么是体系太大跑不动，要么是结果对不上实验吸收谱，要么是根本不知道时间步长该取多少。Octopus 这个开源实空间 TDDFT 代码，恰好是解决这类问题的常用工具之一——它不依赖平面波基组，而是在实空间网格上直接传播含时 Kohn-Sham 方程，天然适合做强场电离、高次谐波、光吸收谱和实时电子动力学。这篇笔记面向的是已经会跑基态 DFT、想把手伸向含时模拟的从业者，也适合做光学响应、分子激发态、纳米结构等方向的人。我会按“概念立住 → 最小可跑 → 参数怎么调 → 坑在哪 → 进阶技巧”的顺序，把 Octopus 做 TDDFT 的落地路径拆开讲，每一步都尽量给到能直接抄的命令和输入片段。

2. Octopus 做 TDDFT 的底层逻辑与最小可跑案例

2.1 为什么是实空间网格，而不是平面波基组

TDDFT 的核心是求解含时 Kohn-Sham 方程：

i ∂ψ_j(r,t)/∂t = [ -∇²/2 + V_KS(r,t) ] ψ_j(r,t)

平面波基组在处理周期性体系时很优雅，但做含时传播时，每个时间步都要做大量 FFT，而且对真空层和孤立分子的描述需要很大的超胞。Octopus 的选择是实空间离散化：把波函数直接存在三维网格上，用有限差分近似动能算符。这样做的好处是：

• 真空区域不需要额外基函数，边界条件用吸收边界或零边界直接处理；
• 时间传播可以用显式算法，每个时间步只做局域乘法；
• 强场下的电离、电子云外溢，实空间网格天然能描述。

代价是网格步长必须足够小，通常要保证每个方向至少 0.2–0.3 Bohr 的分辨率，否则动能算符误差会直接毁掉激发能。常见做法是先用收敛测试确定网格间距，再固定下来做含时计算。

2.2 最小可跑案例：一氧化碳分子的光吸收谱

下面这个例子是我一般用来验证环境是否配好的最小体系：CO 分子，用 LDA 交换关联，先做基态，再做含时传播，最后通过极化率得到吸收谱。输入文件分成两块：基态 inp 和含时 inp。

基态输入inp：

# 基态计算：CO 分子，LDA，实空间网格 CalculationMode = gs UnitsOutput = eV_angstrom %Coordinates "C" | 0.0 | 0.0 | 0.0 "O" | 0.0 | 0.0 | 1.13 % BoxShape = sphere Radius = 6.0 Spacing = 0.25 %Species "C" | species_pseudo | set | standard | lda | 1.0 "O" | species_pseudo | set | standard | lda | 1.0 % TheoryLevel = dft XCFunctional = lda_x + lda_c_pz ConvRelDens = 1e-6 MaximumIter = 200 %Output wfs density %

含时输入inp（覆盖基态文件中的对应块）：

# 含时计算：施加弱电场，计算极化率 CalculationMode = td UnitsOutput = eV_angstrom %Coordinates "C" | 0.0 | 0.0 | 0.0 "O" | 0.0 | 0.0 | 1.13 % BoxShape = sphere Radius = 6.0 Spacing = 0.25 %Species "C" | species_pseudo | set | standard | lda | 1.0 "O" | species_pseudo | set | standard | lda | 1.0 % TheoryLevel = dft XCFunctional = lda_x + lda_c_pz TDPropagator = aetrs TDTimeStep = 0.002 TDMaxSteps = 5000 %TDExternalFields electric_field | 0.01 | 0.0 | 0.0 | envelope_constant | 0.0 | 0.0 % %TDFrequencies 0.0 | 0.0 | 0.0 % %TDOutput multipoles %

运行命令：

# 先跑基态 octopus > gs.log # 再跑含时，注意会读取基态波函数 octopus > td.log

逻辑说明：基态计算生成restart/目录下的波函数和密度；含时计算默认从restart/读取初态。TDExternalFields里给了一个沿 x 方向的常数电场，幅度 0.01 a.u.，属于弱场线性响应区间。TDPropagator = aetrs是近似 enforced time-reversal symmetry 传播子，对大多数体系稳定。TDTimeStep = 0.002对应约 0.05 as，对价电子动力学足够。TDMaxSteps = 5000总时长约 10 fs，频率分辨率约 0.4 eV，对宽谱吸收峰够用，但想看精细振动结构要加长。

参数说明：Spacing是最关键参数，0.25 Å 对 C、O 的赝势够用，但如果用全电子或含 3d 元素，要降到 0.15–0.18。Radius要保证波函数在边界衰减到可忽略，否则吸收谱会出现虚假振荡。TDTimeStep不能太大，一般取 0.001–0.005，太大直接发散。

2.3 从极化率到吸收谱：后处理怎么做

含时跑完后，td.log里会有multipoles输出，记录每个时间步的偶极矩。吸收谱正比于频率乘以极化率虚部。Octopus 自带oct-propagation_spectrum工具：

# 提取吸收谱 oct-propagation_spectrum

它会读取td.general/multipoles，输出cross_section文件。用 gnuplot 或 Python 画出来：

import numpy as np import matplotlib.pyplot as plt data = np.loadtxt('cross_section') plt.plot(data[:,0], data[:,1]) plt.xlabel('Energy (eV)') plt.ylabel('Cross section (arb. units)') plt.xlim(0, 30) plt.show()

逻辑说明：cross_section第一列是能量，第二列是截面。CO 的典型吸收峰在 8–14 eV 之间，如果算出来峰位偏移超过 0.5 eV，优先检查网格间距和交换关联泛函。LDA 通常会低估激发能，这是泛函本身的问题，不是代码 bug。

3. 参数怎么设：网格、时间步、传播子与交换关联的取舍

3.1 网格间距与盒子半径的收敛测试

实空间方法最怕的就是“网格不够细”。我一般会做一组收敛测试：固定盒子半径，逐步减小 Spacing，看总能量和第一激发能的变化。下面是一个批量脚本的骨架：

#!/bin/bash # 收敛测试：不同 Spacing 下的基态总能量 for sp in 0.30 0.25 0.20 0.18 0.15; do mkdir -p sp_$sp cp inp sp_$sp/ sed -i "s/Spacing = .*/Spacing = $sp/" sp_$sp/inp cd sp_$sp octopus > gs.log grep "Total energy" gs.log | tail -1 cd .. done

逻辑说明：每次改 Spacing 后重新跑基态，记录总能量。当相邻两次总能量差小于 1 meV/atom 时，认为收敛。对含时计算，还要看第一激发峰位置是否稳定。经验值：C、H、O、N 用 0.20–0.25 Å；过渡金属用 0.15–0.18 Å；含 f 电子用 0.12–0.15 Å。

盒子半径的测试类似，但更简单：从 4 Å 开始，每次加 1 Å，看总能量和偶极矩是否变化。对弱场吸收谱，半径不够会导致低能端出现虚假峰，因为波函数被边界反射。

3.2 时间步长与传播子的选择

Octopus 支持多种传播子，常用的有：

时间步长的选择有个经验公式：TDTimeStep ≤ 0.1 / E_max，其中 E_max 是体系最高占据态能量的两倍。对价电子体系，0.002–0.005 通常安全。如果跑着跑着能量不守恒，先查时间步长，再查传播子。

# 检查能量守恒：td.log 中 Total energy 应基本恒定 grep "Total energy" td.log | head -20

如果能量漂移超过 1e-4 a.u.，说明时间步长太大或传播子不合适。强场下建议用 etrs，虽然慢一点，但不会翻车。

3.3 交换关联泛函与赝势的匹配

TDDFT 的激发能对泛函很敏感。LDA 通常低估 0.3–0.8 eV，GGA（如 PBE）略好，但仍有系统误差。如果做吸收谱，常见做法是用 LDA 或 PBE 先跑，再用剪刀算符平移。如果做强场电离，泛函影响相对小，因为电离主要由外场决定。

赝势方面，Octopus 自带standard和hgh等系列。注意：赝势的价电子数要和 XCFunctional 匹配。比如用lda赝势就配lda_x + lda_c_pz，不要混用。混用会导致总能量异常，含时传播直接发散。

# 检查赝势文件是否被正确读取 grep "Species" gs.log grep "valence" gs.log

4. 避坑与排查：TDDFT 跑不起来时先看这五条

4.1 含时计算第一步就发散，能量变成 NaN

现象：td.log里 Total energy 在几步内变成 NaN 或极大值。

原因：最常见的是时间步长太大，或者网格间距太粗导致动能算符误差爆炸。其次可能是赝势与泛函不匹配，或者初始波函数没读进去。

解决：先把 TDTimeStep 降到 0.001，Spacing 降到 0.18，重跑。如果还发散，检查restart/目录是否存在，基态是否收敛。基态没收敛就做含时，等于在错误初态上传播，必炸。

4.2 吸收谱出现负值或剧烈振荡

现象：cross_section在低能端出现负截面，或者整条谱线毛刺很多。

原因：总传播时间太短，频率分辨率不够；或者盒子半径太小，波函数被边界反射；或者电场幅度太大，进入非线性区。

解决：把 TDMaxSteps 加到 10000 以上，Radius 加 2 Å，电场幅度降到 0.001。如果还振荡，检查multipoles输出是否等间隔，时间步长是否稳定。

4.3 激发能比实验值低 1 eV 以上

现象：主吸收峰位置明显偏低。

原因：LDA/PBE 的固有误差，不是代码问题。TDDFT 在绝热近似下对电荷转移激发和双激发描述很差。

解决：换用杂化泛函（如 B3LYP）或做剪刀修正。Octopus 支持XCFunctional = b3lyp，但计算量会大很多。如果只是看趋势，LDA 加 0.5 eV 平移也能用。

4.4 并行效率低，跑几天没进展

现象：MPI 并行后加速比远低于线性。

原因：实空间网格的并行域分解在体系较小时通信开销大；或者 FFT 在含时传播中成为瓶颈。

解决：小体系（<50 原子）用 OpenMP 线程并行更划算；大体系用 MPI，但每个域至少要有 10^4 个网格点。检查octopus编译时是否链接了并行 FFTW。

4.5 重启含时计算后结果对不上

现象：中断后从restart/继续跑，谱线形状变了。

原因：含时传播的相位信息存在波函数里，重启时如果时间步编号或外场相位没对齐，会引入跳变。

解决：尽量一次跑完。如果必须重启，用TDMaxSteps分段，但每段都要从同一初态重新激发，不要接着上一段的末态跑。或者用FromScratch = no并确保td.general/multipoles是追加模式。

5. 进阶技巧：用 delta 脉冲和频谱分析快速定位激发态

5.1 delta 脉冲：一次传播拿到全谱

弱场吸收谱的常规做法是加一个常数电场，但更高效的是加 delta 脉冲：在 t=0 时刻给一个瞬时电场，然后自由传播。这样一次计算就能得到整个频率范围的响应。

%TDExternalFields electric_field | 0.001 | 0.0 | 0.0 | envelope_delta | 0.0 | 0.0 %

逻辑说明：envelope_delta在第一个时间步施加电场，之后为零。幅度要小，保证线性响应。传播时间要足够长，频率分辨率 ΔE = 2π / T_max。比如 T_max = 20 fs，ΔE ≈ 0.2 eV。

5.2 用 Python 做峰值拟合与振子强度提取

跑完 delta 脉冲后，multipoles文件里是偶极矩随时间的变化。做傅里叶变换得到极化率，再提取峰位和振子强度：

import numpy as np from scipy.signal import find_peaks # 读取偶极矩 data = np.loadtxt('td.general/multipoles') t = data[:,0] dip = data[:,1] # x 方向偶极矩 # 去趋势，加窗 dip = dip - np.mean(dip) window = np.hanning(len(dip)) dip_fft = np.fft.rfft(dip * window) freq = np.fft.rfftfreq(len(t), d=t[1]-t[0]) # 转成能量，取虚部 energy = freq * 27.2114 # Hartree -> eV strength = -np.imag(dip_fft) # 找峰 peaks, _ = find_peaks(strength, height=0.1*np.max(strength)) for p in peaks: print(f"Peak at {energy[p]:.2f} eV, strength {strength[p]:.3f}")

逻辑说明：加 Hanning 窗减少截断效应；rfft输出复数，虚部对应吸收。峰位就是激发能，峰高正比于振子强度。如果峰位和实验差 0.5 eV 以内，说明泛函和网格都合理。

5.3 一个我常犯的错误：忘了检查盒子边界

早期我做纳米团簇的 TDDFT，吸收谱低能端总有一个莫名其妙的峰。查了两天才发现是盒子半径不够，电子云在边界被反射，形成了驻波。后来养成习惯：每次含时计算前，先看基态密度在边界处的值，如果大于 1e-6，就把 Radius 加 2 Å。这个检查用octopus输出的density文件就能做：

# 用 octopus 自带工具看密度边界值 oct-center-geom < inp

或者直接在gs.log里搜Density at boundary。这个习惯帮我省了很多后悔药。希望帮到你。

本文还有配套的精品资源，点击获取