新建
分子动力学算法及运行实例
满来福@pkuphy
推荐镜像 :DeePMD-kit:v3.0.0a1-2024Q1
推荐机型 :c2_m4_cpu
赞 1
目录
分子动力学算法及运行实例
目录
1. 分子动力学算法
1.1 相空间与正则方程
1.2 相互作用势能，机器学习势函数
1.3 边界条件与系统尺度
1.4 数值分子动力学算法
1.5 系综模拟算法
2. 使用LAMMPS进行分子动力学计算
2.1 软件参数介绍
2.2 运算实例：金属钛在熔点附近的物理性质
2.3 深度势能模块
3. 使用ABACUS进行第一性原理分子动力学计算
3.1 参数设置
3.2 计算性能分析
3.3 结果与展望
4. 参考文献
代码
文本
1 分子动力学算法
分子动力学（MD）是一种数值模拟方法，用于分析原子和分子的物理运动。分子动力学研究原子和分子在一定时间内的相互作用，从而观察系统的动态“演化”。通常，原子和分子的轨迹通过数值求解相互作用粒子系统的运动方程来确定，其中粒子之间的力以及的势能通常使用势函数来计算。
1.1 相空间与正则方程
在经典力学中，无耗散的孤立系统的演化遵从哈密顿正则方程：
其中H是系统的哈密顿量，p和q分别是正则动量和正则坐标。
由正则坐标和正则动量张成的6N维空间称为系统的相空间。系统在某一时刻的状态可以由相空间中的点描述，系统随着时间演化，代表点的坐标也会随时间改变，在向空间中表示为一条轨迹。由正则方程可以得到，相空间中的一点在上一时刻和下一时刻演化的方向是一定的，所以相空间中的轨迹互不交叉（轨迹不与自身交叉，不同轨迹也互不交叉）。
相空间中，体系能量相等的点的集合称为能量曲面。如果系统为孤立体系（能量守恒），则系统的轨迹始终处于同一个能量曲面，不同的能量曲面也互不相交。
根据正则方程可以推出刘维尔定理：，即相空间中一个点邻域内的点的密度不随时间变化。
1.2 相互作用势能，机器学习势函数
原子之间的相互作用势能较为复杂，主要由原子核之间、原子核与电子之间、电子与电子之间的库伦势以及多极势有关。另外，除两体相互作用外，还存在多体相互作用修正项。
1964年，Rahman使用Lennard-Jones势计算了氩原子的分子动力学模型。Lennard-Jones势可以较好地描述稀有气体原子之间的相互作用，但对于其他电子结构更加复杂、非对称的原子与分子，需要更加复杂的势函数，例如EAM多体相互作用势函数等。
原子之间的相互作用势能的起源要追溯到外场中的原子核和电子的波函数计算，因此人们提出了第一性原理分子动力学（AIMD），通过计算电子波函数来获得准确的相互作用势（原子核用经典方法处理）。但这样的计算在每一步都要进行电子自洽迭代计算，其计算量远大于经典分子动力学。之后，还提出了分子动力学的全量子计算，即用量子力学方法处理原子核和电子，计算量更大。这样大的计算量导致第一性原理难以进行长时间的计算，即经典经验势方法“快而不准”，第一性原理“准而不快”。
机器学习势函数的方法部分解决了以上问题。机器学习势函数是一种利用机器学习技术构建的原子间相互作用模型。它通过学习大量第一性原理计算数据，建立起原子结构和能量之间的关系，从而可以高效地预测材料在不同条件下的性质。
使用机器学习得到势函数，首先需要对势函数进行建模，即构造合适的势函数形式。使用深度神经网络对第一性原理分子动力学的计算结果进行拟合，得到拟合的势函数具体形式，通常也称为深度势能。深度势能方法既保留了第一性原理计算的精确性，同时大大简化了计算复杂度，使得原先需要进行自洽迭代计算的势函数可以用较少的参量直接得到。
但目前深度势能方法还存在一些问题，例如模型在不同原子之间的相互作用上的泛化性能不足、自身的对称性不足等。现在普遍使用的方法是对每一种相互作用势能进行分别的建模，并使用数据增强方法在训练过程中引入对称性，从而期望模型隐含对称性。总之，深度势能方法极大地提升了高度准确性下的计算速度，是分子动力学领域的一大突破。
1.3 边界条件与系统尺度
实际的热力学系统中，粒子的数目是很大的，约在1023量级。分子动力学对这样量级的粒子的计算是无能为力的。但处在平衡状态的体系存在一定的平移对称性，即系统的某些参量在较大尺度（远大于原子尺度）上是均匀的。所以，可以选取一个较小的计算范围，采取周期性边界条件来模拟很大的体系。选用周期性边界条件是为了消除表面效应（例如处于固体/液体表面的分子的行为与处于内部的分子不同）。
然而，用小体系加上周期性边界条件的计算结果与大体系的计算结果存在一些偏差。由于周期性边界条件会引入长程关联，所以有些物理量（例如径向分布函数和均方位移）会与实际情况出现偏差。
总体上，周期性边界条件的方法还是为较大体系的模拟提供了有效的方法，可以用102~104量级的粒子基本模拟1023量级的粒子的行为。
1.4 数值分子动力学算法
根据正则方程，分子动力学的计算的本质是一系列常微分方程的前向求解。在分子动力学领域，最普遍使用的算法是Verlet算法和Velocity-Verlet算法。
Verlet算法的核心是使用当前时刻和上一时刻的坐标，以及当前时刻的力（即势函数的梯度）得到下一时刻的坐标：
Velocity-Verlet算法则是增加了速度信息作为输出量：
可以证明，每一步粒子坐标r的误差在量级，粒子速度v的误差在量级。由于粒子坐标的误差会随时间积累，所以长时间计算仍然会导致很大的误差。根据微分方程的李雅普诺夫不稳定性，不同的系统即使初始条件相差很小，随着时间的演化，其差异也会随时间指数上升。有人认为分子动力学的数值计算结果可能代表某一个系统的演化（即shadow orbitals假说）。总之，分子动力学模拟的是多体系统的统计行为，一般情况关注的点多为系统的统计量，而统计量对此类微扰是不敏感的。
1.5 系综模拟算法
统计物理中，系综代表一定条件下一个体系的大量可能微观状态的集合。分子动力学计算的结果实质是对系统处于平衡态（或附近）的状态进行采样，并认为其频率分布等同于实际系综中的概率分布。
为了使系统约束在指定的外界条件（例如温度、压强、能量）下，需要使用不同的系综模拟算法。下面首先介绍几种系综：
- 微正则系综（NVE系综）：系统的粒子数不变，体积不变，总能量不变（即孤立体系）；
- 正则系综（NVT系综）：系统的粒子数不变，体积不变，温度不变（相当于系统与大热源接触，与大热源交换热量达到平衡）；
- 等温等压系综（NPT系综）：系统的粒子数不变，压强不变，温度不变（相当于系统同时与大热源和恒压热源接触，与大热源交换热量，与恒压热源达到受力平衡）；
- 巨正则系综（μVT系综）：系统的粒子数可变，但化学势（系统增加单个粒子增加的能量）恒定，温度不变（相当于系统与大热源接触，同时可以自由交换粒子）。
模拟这些系综的算法是在系统的拉氏量中增加与守恒量相关的自由度，粒子的运动状态与这些自由度也存在一定的关系。每次演化时同时演化这些附加自由度，通过特定的方法来实现统计参量的恒定。以正则系综为例，体积的固定是通过边界条件不变，温度的固定则存在Langevin算法、Nosé-Hoover算法等方法。Nosé-Hoover算法是典型的通过引入额外自由度，通过衰减这些自由度使体系温度（速率分布函数）稳定。Langevin算法基于布朗运动的随机力理论，通过引入与温度相关的随机力使粒子的速率分布稳定。
代码
文本
2 使用LAMMPS进行分子动力学计算
现在，市面上已经存在多种高效的分子动力学数值模拟软件。下面将以LAMMPS软件为例，介绍分子动力学的计算过程。
2.1 软件参数介绍
LAMMPS 是一款经典的分子动力学（MD）软件，用于模拟液体、固体或气体状态的粒子集合。它可以使用各种原子间势函数和边界条件来模拟原子、聚合物、生物、固体、颗粒或宏观系统。它可以模拟 2D 或 3D 系统，其大小范围从只有几个粒子到数十亿个粒子不等。它支持多种力场和边界条件，并具有高度的并行计算能力。LAMMPS 还易于修改和扩展，可以满足不同的研究需求。
LAMMPS是一款功能强大的软件，可以自定义多种边界条件、势函数、系综、约束、体系结构和计算方式。
LAMMPS至少需要一个输入文件，如果需要导入势函数，还需要势函数文件。下面将简要介绍LAMMPS的输入文件。
2.2 运算实例：金属钛在熔点附近的物理性质
下面给出了金属钛的分子动力学计算输入文件：
# initialize simulation settings
units metal
boundary p p p
atom_style atomic
# define the simulation cell
lattice bcc 3.32
region box block 0 5 0 5 0 5
create_box 1 box
create_atoms 1 box
group Ti type 1
mass 1 47.867
# set force field
pair_style deepmd Ti.pb
pair_coeff * *
# npt simulation
velocity all create 1600 23456789
fix 1 all npt temp 1600 1600 0.5 iso 1.0 1.0 0.5
timestep 0.001
# rdf calculation
compute myRDF all rdf 50
fix 2 all ave/time 100 1 100 c_myRDF[*] file tmp1.rdf1 mode vector
# msd calculation
compute myMSD all msd
fix 3 all ave/time 100 1 100 c_myMSD[*] file Ti.msd1
# output
thermo_style custom step temp pe ke etotal press lx ly lz vol
thermo 100
dump 1 all custom 100 Ti1.dump id type x y z
log log1.lammps
run 5000
#号后的内容是注释。
units代表计算采用的单位制，metal单位使用g/mol作为质量单位，Angstrom作为长度单位，ps作为时间单位，开尔文作为温度单位等等。
boundary代表边界条件，p表示周期性边界条件。
lattice 代表原子的初始位置使用晶格放置，bcc代表晶格种类，数字为晶格常数。后面的几行代表在每一个格点创建一个钛原子，region后的参数代表创建5×5×5的网格（即一个超胞含有125个单胞，250个原子）。
pair_style表示采用的势函数。这里导入了deepmd模块，并指定了导入的势函数文件。
pair_coeff代表对势的系数，这里为缺省值。
velocity参数控制系统的速度初始化。其格式为：
velocity group-ID style args keyword value ...
其中，style可取create、set、scale、ramp、zero等。create表示按速率分布创建初始速度分布，使用了温度1600K、随机种子23456789。
timestep表示时间步长。
compute表示输出的计算参数。
run表示计算指定步数。
fix表示指定某些参数，例如npt表示指定NPT系综，ave/time表示将时间平均的计算结果输出至指定文件。
下面是一个模型的运行结果（由于计算成本限制，这里略去了其他温度下的运行，下面直接展示在不同温度下计算的几次数据）：
代码
文本
[1]
%%bash
pwd
cd personal/lammps_learn/lab3/
lmp -i Ti.in
/ LAMMPS (2 Aug 2023 - Update 3) OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (src/comm.cpp:98) using 1 OpenMP thread(s) per MPI task DeePMD-kit: Successfully load libcudart.so.12 DeePMD-kit WARNING: Environmental variable DP_INTRA_OP_PARALLELISM_THREADS is not set. Tune DP_INTRA_OP_PARALLELISM_THREADS for the best performance. See https://deepmd.rtfd.io/parallelism/ for more information. DeePMD-kit WARNING: Environmental variable DP_INTER_OP_PARALLELISM_THREADS is not set. Tune DP_INTER_OP_PARALLELISM_THREADS for the best performance. See https://deepmd.rtfd.io/parallelism/ for more information. DeePMD-kit WARNING: Environmental variable OMP_NUM_THREADS is not set. Tune OMP_NUM_THREADS for the best performance. See https://deepmd.rtfd.io/parallelism/ for more information. Loaded 1 plugins from /opt/deepmd-kit/dp/lib/deepmd_lmp Lattice spacing in x,y,z = 3.32 3.32 3.32 Created orthogonal box = (0 0 0) to (16.6 16.6 16.6) 1 by 1 by 1 MPI processor grid Created 250 atoms using lattice units in orthogonal box = (0 0 0) to (16.6 16.6 16.6) create_atoms CPU = 0.000 seconds 250 atoms in group Ti Summary of lammps deepmd module ... >>> Info of deepmd-kit: installed to: /opt/deepmd-kit/dp source: v3.0.0a0-83-gd23cf3e6 source branch: mapping source commit: d23cf3e6 source commit at: 2024-05-07 06:32:03 +0000 support model ver.: 1.1 build variant: cuda build with tf inc: /opt/mamba/lib/python3.10/site-packages/tensorflow/include;/opt/mamba/lib/python3.10/site-packages/tensorflow/include build with tf lib: /opt/mamba/lib/python3.10/site-packages/tensorflow/libtensorflow_cc.so.2 build with pt lib: torch;torch_library;/opt/libtorch/lib/libc10.so;/opt/libtorch/lib/libkineto.a;/usr/local/cuda/targets/x86_64-linux/lib/stubs/libcuda.so;/usr/local/cuda-12.1/targets/x86_64-linux/lib/libnvrtc.so;/usr/local/cuda/lib64/libnvToolsExt.so;/usr/local/cuda/lib64/libcudart.so;/opt/libtorch/lib/libc10_cuda.so set tf intra_op_parallelism_threads: 0 set tf inter_op_parallelism_threads: 0 >>> Info of lammps module: use deepmd-kit at: /opt/deepmd-kit/dp source: v3.0.0a0-83-gd23cf3e6 source branch: mapping source commit: d23cf3e6 source commit at: 2024-05-07 06:32:03 +0000 build float prec: double build with tf inc: /opt/mamba/lib/python3.10/site-packages/tensorflow/include;/opt/mamba/lib/python3.10/site-packages/tensorflow/include build with tf lib: /opt/mamba/lib/python3.10/site-packages/tensorflow/libtensorflow_cc.so.2 DeePMD-kit WARNING: Environmental variable DP_INTRA_OP_PARALLELISM_THREADS is not set. Tune DP_INTRA_OP_PARALLELISM_THREADS for the best performance. See https://deepmd.rtfd.io/parallelism/ for more information. DeePMD-kit WARNING: Environmental variable DP_INTER_OP_PARALLELISM_THREADS is not set. Tune DP_INTER_OP_PARALLELISM_THREADS for the best performance. See https://deepmd.rtfd.io/parallelism/ for more information. DeePMD-kit WARNING: Environmental variable OMP_NUM_THREADS is not set. Tune OMP_NUM_THREADS for the best performance. See https://deepmd.rtfd.io/parallelism/ for more information. 2024-06-03 08:44:06.107572: I tensorflow/core/platform/cpu_feature_guard.cc:210] This TensorFlow binary is optimized to use available CPU instructions in performance-critical operations. To enable the following instructions: AVX2 AVX512F FMA, in other operations, rebuild TensorFlow with the appropriate compiler flags. 2024-06-03 08:44:06.175972: I tensorflow/compiler/mlir/mlir_graph_optimization_pass.cc:388] MLIR V1 optimization pass is not enabled INVALID_ARGUMENT: Tensor spin_attr/ntypes_spin:0, specified in either feed_devices or fetch_devices was not found in the Graph >>> Info of model(s): using 1 model(s): Ti.pb rcut in model: 9 ntypes in model: 1 CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE Your simulation uses code contributions which should be cited: - USER-DEEPMD package: The log file lists these citations in BibTeX format. CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule Neighbor list info ... update: every = 1 steps, delay = 0 steps, check = yes max neighbors/atom: 2000, page size: 100000 master list distance cutoff = 11 ghost atom cutoff = 11 binsize = 5.5, bins = 4 4 4 2 neighbor lists, perpetual/occasional/extra = 1 1 0 (1) pair deepmd, perpetual attributes: full, newton on pair build: full/bin/atomonly stencil: full/bin/3d bin: standard (2) compute rdf, occasional, half/full from (1) attributes: half, newton on pair build: halffull/newton stencil: none bin: none Setting up Verlet run ... Unit style : metal Current step : 0 Time step : 0.001 Per MPI rank memory allocation (min/avg/max) = 4.718 | 4.718 | 4.718 Mbytes Step Temp PotEng KinEng TotEng Press Lx Ly Lz Volume 0 1600 -1925.8177 51.497242 -1874.3204 -47840.185 16.6 16.6 16.6 4574.296 100 1063.6922 -1908.5127 34.235758 -1874.277 18345.006 16.350663 16.350663 16.350663 4371.2544 200 1011.7507 -1905.6211 32.563981 -1873.0571 -14552.325 16.544227 16.544227 16.544227 4528.344 300 1140.0403 -1908.0516 36.693083 -1871.3585 13003.954 16.398843 16.398843 16.398843 4410.0109 400 1099.7779 -1905.2529 35.397207 -1869.8557 -10049.167 16.526862 16.526862 16.526862 4514.0999 500 1170.0331 -1906.1043 37.658423 -1868.4459 6645.7904 16.422319 16.422319 16.422319 4428.9771 600 1212.635 -1905.8241 39.029598 -1866.7945 -6752.4715 16.517202 16.517202 16.517202 4506.1893 700 1278.7063 -1906.1547 41.156156 -1864.9986 4214.3347 16.442696 16.442696 16.442696 4445.4847 800 1220.2555 -1902.2817 39.27487 -1863.0069 -4356.9199 16.507664 16.507664 16.507664 4498.3878 900 1264.3389 -1901.6197 40.693729 -1860.926 4341.0658 16.462632 16.462632 16.462632 4461.6739 1000 1293.4894 -1900.4635 41.631962 -1858.8315 -5748.1519 16.518516 16.518516 16.518516 4507.2651 1100 1354.1912 -1900.3481 43.585695 -1856.7624 5751.0329 16.476033 16.476033 16.476033 4472.5785 1200 1436.5152 -1900.93 46.235358 -1854.6946 -6654.0873 16.538766 16.538766 16.538766 4523.8619 1300 1289.1879 -1894.0442 41.493513 -1852.5507 4837.1834 16.467561 16.467561 16.467561 4465.6821 1400 1429.0187 -1896.4979 45.994075 -1850.5038 -5583.4138 16.524333 16.524333 16.524333 4512.028 1500 1430.7793 -1894.4961 46.050743 -1848.4453 2332.335 16.491565 16.491565 16.491565 4485.2396 1600 1501.9133 -1894.8198 48.340245 -1846.4795 -7151.0868 16.551141 16.551141 16.551141 4534.0236 1700 1473.5832 -1892.0094 47.42842 -1844.581 8281.299 16.486556 16.486556 16.486556 4481.1533 1800 1552.801 -1892.697 49.978106 -1842.7189 -12187.2 16.601809 16.601809 16.601809 4575.7915 1900 1659.3419 -1894.3865 53.407207 -1840.9793 14040.893 16.458735 16.458735 16.458735 4458.5058 2000 1590.3287 -1890.5876 51.185962 -1839.4016 -19398.397 16.653891 16.653891 16.653891 4618.9917 2100 1665.5373 -1891.4451 53.60661 -1837.8385 23679.361 16.420589 16.420589 16.420589 4427.578 2200 1638.912 -1889.3793 52.749655 -1836.6297 -21783.457 16.677259 16.677259 16.677259 4638.462 2300 1650.6243 -1888.8188 53.126623 -1835.6921 18643.552 16.450576 16.450576 16.450576 4451.8788 2400 1682.8108 -1889.1283 54.162571 -1834.9657 -15155.224 16.644216 16.644216 16.644216 4610.9458 2500 1727.7936 -1889.9476 55.610377 -1834.3372 7942.7968 16.495366 16.495366 16.495366 4488.3416 2600 1753.7916 -1890.3485 56.447144 -1833.9014 -9669.0635 16.613304 16.613304 16.613304 4585.303 2700 1715.4343 -1888.8911 55.212584 -1833.6785 7922.0966 16.509788 16.509788 16.509788 4500.1243 2800 1811.2499 -1891.8234 58.296483 -1833.5269 -9088.3444 16.587273 16.587273 16.587273 4563.7831 2900 1625.4631 -1885.9571 52.316792 -1833.6403 8325.2473 16.507821 16.507821 16.507821 4498.5162 3000 1730.4671 -1889.6268 55.696427 -1833.9303 -9684.1669 16.59662 16.59662 16.59662 4571.5022 3100 1750.2022 -1890.8305 56.331615 -1834.4989 4599.7829 16.517584 16.517584 16.517584 4506.5021 3200 1644.5681 -1888.3593 52.931701 -1835.4276 -5018.0797 16.587395 16.587395 16.587395 4563.8836 3300 1672.6409 -1890.2409 53.835245 -1836.4057 6468.2836 16.51259 16.51259 16.51259 4502.4154 3400 1629.4813 -1890.1389 52.44612 -1837.6927 -5641.8422 16.590415 16.590415 16.590415 4566.3771 3500 1686.5105 -1893.1755 54.281651 -1838.8939 3737.2917 16.522858 16.522858 16.522858 4510.82 3600 1570.1517 -1890.661 50.536551 -1840.1245 -1682.3715 16.561323 16.561323 16.561323 4542.3972 3700 1541.2649 -1890.7741 49.606807 -1841.1672 -1807.2537 16.541352 16.541352 16.541352 4525.9843 3800 1512.7999 -1890.8019 48.69064 -1842.1112 1845.2596 16.52465 16.52465 16.52465 4512.2881 3900 1658.1374 -1896.2045 53.368438 -1842.836 -3100.4726 16.549638 16.549638 16.549638 4532.7891 4000 1493.7646 -1891.6414 48.077973 -1843.5634 4140.4737 16.510995 16.510995 16.510995 4501.1114 4100 1443.3568 -1890.6106 46.455558 -1844.155 -4462.0918 16.539584 16.539584 16.539584 4524.5329 4200 1503.9992 -1892.829 48.407382 -1844.4216 7184.2847 16.480446 16.480446 16.480446 4476.173 4300 1521.684 -1893.6165 48.976581 -1844.6399 -6462.8769 16.569085 16.569085 16.569085 4548.7868 4400 1541.5616 -1894.2182 49.616357 -1844.6018 3207.0496 16.509348 16.509348 16.509348 4499.7641 4500 1490.2613 -1892.3508 47.965217 -1844.3856 -6008.0499 16.566897 16.566897 16.566897 4546.9847 4600 1600.63 -1895.4933 51.517518 -1843.9758 2674.0593 16.505993 16.505993 16.505993 4497.0218 4700 1610.6413 -1895.2559 51.839741 -1843.4162 -7562.5733 16.579031 16.579031 16.579031 4556.9831 4800 1536.4351 -1892.1049 49.451355 -1842.6536 8383.8232 16.488978 16.488978 16.488978 4483.129 4900 1622.0844 -1893.987 52.208044 -1841.7789 -10401.309 16.576107 16.576107 16.576107 4554.5727 5000 1619.922 -1892.9732 52.138447 -1840.8347 9703.5474 16.479574 16.479574 16.479574 4475.4631 Loop time of 780.153 on 1 procs for 5000 steps with 250 atoms Performance: 0.554 ns/day, 43.342 hours/ns, 6.409 timesteps/s, 1.602 katom-step/s 187.3% CPU use with 1 MPI tasks x 1 OpenMP threads MPI task timing breakdown: Section | min time | avg time | max time |%varavg| %total --------------------------------------------------------------- Pair | 779.54 | 779.54 | 779.54 | 0.0 | 99.92 Neigh | 0.20554 | 0.20554 | 0.20554 | 0.0 | 0.03 Comm | 0.13465 | 0.13465 | 0.13465 | 0.0 | 0.02 Output | 0.018743 | 0.018743 | 0.018743 | 0.0 | 0.00 Modify | 0.21758 | 0.21758 | 0.21758 | 0.0 | 0.03 Other | | 0.04024 | | | 0.01 Nlocal: 250 ave 250 max 250 min Histogram: 1 0 0 0 0 0 0 0 0 0 Nghost: 2815 ave 2815 max 2815 min Histogram: 1 0 0 0 0 0 0 0 0 0 Neighs: 38826 ave 38826 max 38826 min Histogram: 1 0 0 0 0 0 0 0 0 0 FullNghs: 77652 ave 77652 max 77652 min Histogram: 1 0 0 0 0 0 0 0 0 0 Total # of neighbors = 77652 Ave neighs/atom = 310.608 Neighbor list builds = 86 Dangerous builds = 0 Total wall time: 0:13:06
代码
文本
可以通过结果计算径向分布函数（RDF）和均方位移（MSD）：
代码
文本
[8]
import numpy as np
import matplotlib.pyplot as plt
import os
os.chdir('/personal/lammps_learn/lab3/')
msd_files = ['Ti.msd1','Ti2.msd2','Ti3.msd3','Ti4.msd4','Ti5.msd5']
T_list = [1600, 1800, 2000, 2200, 2400]
for filename, T in zip(msd_files, T_list):
data = np.loadtxt(filename, skiprows=2)
time = data[:, 0]
time -= time[0]
msd1 = data[:, 3]
if T>2000:
l='dashed'
else:
l='solid'
plt.plot(time / 1000, msd1, label=f'T={T}k', linestyle=l)
plt.xlabel('time(ps)')
plt.ylabel('$MSD(Å^2)$')
plt.title('msd')
plt.legend()
plt.show()
nbins = 50 # define the number of bins in the RDF
rdf_files = ['tmp1.rdf1','tmp2.rdf2','tmp3.rdf3','tmp4.rdf4','tmp5.rdf5']
for filename, T in zip(rdf_files, T_list):
with open(filename, "r") as f:
lines = f.readlines()
lines = lines[3:]
data = np.zeros((nbins, 3))
count = 0
fact = 0
# 读取文件
for line in lines:
nums = line.split()
if len(nums) == 4 and count > 25:
for i in range(1, 4):
data[int(nums[0]) - 1, i - 1] += float(nums[i])
if len(nums) == 2:
count += 1
# 去掉前20组数据，因为还没有平衡
if count > 20:
fact += 1
ave_rdf = data / fact # 计算平均
if T>2000:
l='dashed'
else:
l='solid'
plt.plot(ave_rdf[:, 0], ave_rdf[:, 1], label=f'T={T}k', linestyle=l)
plt.title('RDF')
plt.xlabel('r/Å')
plt.ylabel('g(r)')
plt.legend()
plt.show()
代码
文本
从运行结果可以看出，在2000K到2200K之间，径向分布函数和均方位移的行为都发生了改变。 在2000K及以下，径向分布函数的几个峰较为尖锐，均方位移没有显著增加，反映了固体中晶格对称性以及原子在原地振动的行为；在2200K及以上，径向分布函数的峰变得平缓，且每个峰右移，反映了液体中短程有序、扩散的特性，且相变时发生密度变化。
代码
文本
3 使用ABACUS进行第一性原理分子动力学计算
下面使用ABACUS中的第一性原理分子动力学计算重复以上步骤：
3.1 参数设置
参数设置如下：
INPUT_PARAMETERS
#Parameters (1.General)
suffix Ti_nvt
calculation md
nbands 100
symmetry 0
pseudo_dir ./
orbital_dir ./
#Parameters (2.Iteration)
ecutwfc 30
scf_thr 1e-4
scf_nmax 100
#Parameters (3.Basis)
basis_type lcao
ks_solver genelpa
gamma_only 1
#Parameters (4.Smearing)
smearing_method gaussian
smearing_sigma 0.001
#Parameters (5.Mixing)
mixing_type pulay
mixing_beta 0.3
chg_extrap second-order
#Parameters (6.MD)
md_type nvt
md_nstep 20
md_dt 1
md_tfirst 2000
md_tfreq 1
md_tchain 1
md_dumpfreq 1
out_stru 1
为了节省计算成本，这里只计算20步，系统大小采用2×2×2个单胞，共16个原子。
代码
文本
[1]
%%bash
cd personal/abacus_learn/lab11_md
mpirun -n 2 abacus
/
ABACUS v3.6.0
Atomic-orbital Based Ab-initio Computation at UStc
Website: http://abacus.ustc.edu.cn/
Documentation: https://abacus.deepmodeling.com/
Repository: https://github.com/abacusmodeling/abacus-develop
https://github.com/deepmodeling/abacus-develop
Commit: 5b263cbff (Wed Mar 27 18:42:03 2024 +0800)
Mon Jun 3 08:12:42 2024
MAKE THE DIR : OUT.Ti_nvt/
MAKE THE STRU DIR : OUT.Ti_nvt/STRU/
RUNNING WITH DEVICE : CPU / Intel(R) Xeon(R) Platinum
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Warning: the number of valence electrons in pseudopotential > 4 for Ti: [Ar] 3d2 4s2
Pseudopotentials with additional electrons can yield (more) accurate outcomes, but may be less efficient.
If you're confident that your chosen pseudopotential is appropriate, you can safely ignore this warning.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
UNIFORM GRID DIM : 45 * 45 * 45
UNIFORM GRID DIM(BIG) : 15 * 15 * 15
DONE(0.596418 SEC) : SETUP UNITCELL
DONE(0.616604 SEC) : INIT K-POINTS
---------------------------------------------------------
Molecular Dynamics simulations
---------------------------------------------------------
ENSEMBLE : NVT mode: nhc
Time interval(fs) : 1
---------------------------------------------------------
SPIN KPOINTS PROCESSORS NBASE
1 Gamma 2 432
---------------------------------------------------------
Use Systematically Improvable Atomic bases
---------------------------------------------------------
ELEMENT ORBITALS NBASE NATOM XC
Ti 4s2p2d1f-8au 27 16
---------------------------------------------------------
Initial plane wave basis and FFT box
---------------------------------------------------------
DONE(0.641991 SEC) : INIT PLANEWAVE
----------------------------------- INIT VEL ---------------------------------------
RANDOM VEL ACCORDING TO INITIAL TEMPERATURE: 2000 K
------------------------------------- DONE -----------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 0
-------------------------------------------
START CHARGE : atomic
DONE(1.13234 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511205e+04 0.000000e+00 2.076e-02 1.464e+01
GE2 -2.511155e+04 4.997518e-01 1.239e-02 3.168e+00
GE3 -2.511157e+04 -1.568636e-02 1.947e-02 3.321e+00
GE4 -2.511171e+04 -1.407906e-01 1.286e-02 3.166e+00
GE5 -2.511171e+04 1.933734e-03 4.578e-03 3.222e+00
GE6 -2.511172e+04 -1.485039e-02 6.096e-03 3.190e+00
GE7 -2.511173e+04 -9.078252e-03 5.848e-03 3.180e+00
GE8 -2.511171e+04 2.143960e-02 6.751e-03 3.170e+00
GE9 -2.511173e+04 -1.511172e-02 2.650e-03 3.140e+00
GE10 -2.511172e+04 1.385025e-03 2.169e-03 3.150e+00
GE11 -2.511172e+04 -8.968681e-04 2.612e-03 3.137e+00
GE12 -2.511172e+04 -7.337556e-05 2.443e-03 3.125e+00
GE13 -2.511172e+04 5.810796e-04 2.771e-03 3.130e+00
GE14 -2.511172e+04 -1.951953e-04 1.142e-03 3.123e+00
GE15 -2.511172e+04 3.321681e-04 1.597e-03 3.163e+00
GE16 -2.511172e+04 -1.880340e-04 3.247e-04 3.153e+00
GE17 -2.511172e+04 4.633817e-05 8.226e-05 3.123e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.850e-01 2.000e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 1
-------------------------------------------
DONE(1.041e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511173e+04 0.000000e+00 3.339e-03 1.450e+01
GE2 -2.511172e+04 5.079691e-03 8.634e-03 3.104e+00
GE3 -2.511173e+04 -3.654656e-03 5.188e-03 3.097e+00
GE4 -2.511173e+04 -7.116994e-03 5.644e-03 3.099e+00
GE5 -2.511174e+04 -2.209829e-03 4.445e-03 3.094e+00
GE6 -2.511174e+04 -3.189736e-03 1.814e-03 3.119e+00
GE7 -2.511174e+04 -7.910661e-04 1.514e-03 3.090e+00
GE8 -2.511174e+04 -4.932853e-04 1.024e-03 3.140e+00
GE9 -2.511174e+04 -2.539302e-04 8.153e-04 3.093e+00
GE10 -2.511174e+04 1.587756e-04 1.060e-03 3.108e+00
GE11 -2.511174e+04 -2.139977e-04 2.678e-04 3.105e+00
GE12 -2.511174e+04 -1.112614e-05 1.730e-04 3.106e+00
GE13 -2.511174e+04 2.988606e-06 2.103e-04 3.096e+00
GE14 -2.511174e+04 2.559317e-07 2.625e-04 3.102e+00
GE15 -2.511174e+04 4.486925e-07 2.322e-04 3.101e+00
GE16 -2.511174e+04 2.648602e-05 8.740e-05 3.094e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.857e-01 2.005e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 2
-------------------------------------------
DONE(2.023e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511174e+04 0.000000e+00 3.377e-03 1.451e+01
GE2 -2.511172e+04 1.820448e-02 1.378e-02 3.144e+00
GE3 -2.511175e+04 -3.309824e-02 5.713e-03 3.144e+00
GE4 -2.511176e+04 -1.065960e-02 2.313e-03 3.102e+00
GE5 -2.511176e+04 -1.416638e-03 1.363e-03 3.133e+00
GE6 -2.511176e+04 -3.740122e-04 1.254e-03 3.205e+00
GE7 -2.511176e+04 1.388109e-05 9.602e-04 3.149e+00
GE8 -2.511176e+04 1.130837e-04 7.200e-04 3.175e+00
GE9 -2.511176e+04 -2.054230e-04 5.348e-04 3.147e+00
GE10 -2.511176e+04 -3.975386e-05 3.124e-04 3.151e+00
GE11 -2.511176e+04 1.288784e-07 3.384e-04 3.134e+00
GE12 -2.511176e+04 9.382963e-06 3.407e-04 3.137e+00
GE13 -2.511176e+04 -2.614421e-05 2.401e-04 3.133e+00
GE14 -2.511176e+04 -2.821253e-06 1.939e-04 3.208e+00
GE15 -2.511176e+04 -5.429714e-06 1.025e-04 3.236e+00
GE16 -2.511176e+04 -2.581130e-06 1.437e-04 3.182e+00
GE17 -2.511176e+04 3.193679e-05 6.599e-05 3.163e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.858e-01 2.005e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 3
-------------------------------------------
DONE(3.051e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511170e+04 0.000000e+00 5.389e-03 1.454e+01
GE2 -2.511168e+04 1.611480e-02 1.679e-02 3.112e+00
GE3 -2.511177e+04 -8.882308e-02 5.027e-03 3.097e+00
GE4 -2.511178e+04 -1.185479e-02 2.844e-03 3.101e+00
GE5 -2.511178e+04 -2.123400e-03 1.704e-03 3.100e+00
GE6 -2.511179e+04 -4.740230e-04 9.556e-04 3.130e+00
GE7 -2.511179e+04 2.812863e-04 8.556e-04 3.198e+00
GE8 -2.511179e+04 -1.600149e-05 7.336e-04 3.170e+00
GE9 -2.511179e+04 -9.442984e-05 4.647e-04 3.160e+00
GE10 -2.511179e+04 -2.416712e-05 4.698e-04 3.134e+00
GE11 -2.511179e+04 -1.304998e-05 4.342e-04 3.169e+00
GE12 -2.511179e+04 -1.475337e-05 2.545e-04 3.134e+00
GE13 -2.511179e+04 -1.930618e-05 2.095e-04 3.128e+00
GE14 -2.511179e+04 -7.001539e-06 1.585e-04 3.148e+00
GE15 -2.511179e+04 9.072028e-06 6.008e-05 3.116e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.844e-01 1.996e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 4
-------------------------------------------
DONE(4.017e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511178e+04 0.000000e+00 2.870e-03 1.467e+01
GE2 -2.511177e+04 8.363093e-03 9.329e-03 3.185e+00
GE3 -2.511180e+04 -3.503134e-02 2.056e-03 3.181e+00
GE4 -2.511180e+04 -2.022496e-03 8.968e-04 3.149e+00
GE5 -2.511180e+04 -2.503391e-04 7.385e-04 3.131e+00
GE6 -2.511180e+04 -8.753110e-05 4.503e-04 3.166e+00
GE7 -2.511180e+04 6.164287e-06 4.248e-04 3.191e+00
GE8 -2.511180e+04 3.112301e-05 3.601e-04 3.194e+00
GE9 -2.511180e+04 3.546772e-05 2.668e-04 3.271e+00
GE10 -2.511180e+04 -1.927902e-05 1.173e-04 3.183e+00
GE11 -2.511180e+04 1.395478e-05 5.008e-05 3.146e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.841e-01 1.994e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 5
-------------------------------------------
DONE(4.856e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511179e+04 0.000000e+00 3.030e-03 1.442e+01
GE2 -2.511177e+04 1.990292e-02 1.225e-02 3.117e+00
GE3 -2.511182e+04 -4.737823e-02 1.374e-03 3.122e+00
GE4 -2.511182e+04 -8.832862e-04 7.671e-04 3.119e+00
GE5 -2.511182e+04 -2.363218e-04 5.307e-04 3.116e+00
GE7 -2.511182e+04 1.296443e-05 4.040e-04 3.133e+00
GE8 -2.511182e+04 1.228436e-05 3.211e-04 3.118e+00
GE9 -2.511182e+04 1.530486e-05 2.121e-04 3.142e+00
GE10 -2.511182e+04 -1.744407e-06 3.581e-05 3.113e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.851e-01 2.001e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 6
-------------------------------------------
DONE(5.647e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511177e+04 0.000000e+00 4.818e-03 1.421e+01
GE2 -2.511175e+04 1.771166e-02 1.186e-02 3.167e+00
GE3 -2.511183e+04 -7.661979e-02 4.671e-03 3.144e+00
GE4 -2.511183e+04 -7.090926e-03 2.195e-03 3.128e+00
GE5 -2.511183e+04 -2.136220e-03 1.860e-03 3.143e+00
GE6 -2.511183e+04 -1.712357e-04 1.422e-03 3.104e+00
GE7 -2.511183e+04 3.540956e-04 1.309e-03 3.177e+00
GE8 -2.511183e+04 3.030190e-07 6.986e-04 3.150e+00
GE9 -2.511183e+04 1.090617e-04 4.732e-04 3.182e+00
GE10 -2.511183e+04 -6.138468e-05 1.843e-04 3.160e+00
GE11 -2.511183e+04 4.353744e-05 6.967e-05 3.105e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.856e-01 2.004e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 7
-------------------------------------------
DONE(6.457e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511181e+04 0.000000e+00 3.262e-03 1.412e+01
GE2 -2.511180e+04 7.114662e-03 1.072e-02 3.247e+00
GE3 -2.511184e+04 -4.032815e-02 3.565e-03 3.214e+00
GE4 -2.511184e+04 -4.524729e-03 9.541e-04 3.158e+00
GE5 -2.511184e+04 -4.410958e-04 5.009e-04 3.133e+00
GE6 -2.511184e+04 -2.046807e-05 3.337e-04 3.140e+00
GE7 -2.511184e+04 7.934200e-06 3.716e-04 3.231e+00
GE8 -2.511184e+04 2.990552e-05 2.832e-04 3.196e+00
GE9 -2.511184e+04 -1.286619e-06 1.306e-04 3.283e+00
GE10 -2.511184e+04 1.295670e-05 4.591e-05 3.296e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.850e-01 2.000e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 8
-------------------------------------------
DONE(7.247e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511180e+04 0.000000e+00 3.656e-03 1.435e+01
GE2 -2.511180e+04 4.135093e-03 1.139e-02 3.223e+00
GE3 -2.511185e+04 -4.783062e-02 3.671e-03 3.188e+00
GE4 -2.511185e+04 -5.259838e-03 1.110e-03 3.215e+00
GE5 -2.511185e+04 -6.301624e-04 7.053e-04 3.200e+00
GE6 -2.511185e+04 -3.386755e-05 4.029e-04 3.204e+00
GE7 -2.511185e+04 -9.534019e-06 4.707e-04 3.226e+00
GE8 -2.511185e+04 4.431655e-05 3.927e-04 3.181e+00
GE9 -2.511185e+04 4.671287e-06 2.016e-04 3.204e+00
GE10 -2.511185e+04 2.717197e-05 5.879e-05 3.172e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.844e-01 1.995e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 9
-------------------------------------------
DONE(8.025e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511183e+04 0.000000e+00 2.142e-03 1.432e+01
GE2 -2.511184e+04 -1.240716e-03 5.278e-03 3.201e+00
GE3 -2.511185e+04 -1.596001e-02 2.070e-03 3.211e+00
GE4 -2.511185e+04 -1.857701e-03 4.430e-04 3.162e+00
GE5 -2.511185e+04 -8.385508e-05 2.732e-04 3.141e+00
GE6 -2.511185e+04 -1.266796e-05 1.627e-04 3.148e+00
GE7 -2.511185e+04 7.651067e-06 1.491e-04 3.170e+00
GE8 -2.511185e+04 3.794635e-06 1.236e-04 3.223e+00
GE9 -2.511185e+04 2.081560e-05 5.502e-05 3.219e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.846e-01 1.997e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 10
-------------------------------------------
DONE(8.773e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511178e+04 0.000000e+00 4.707e-03 1.416e+01
GE2 -2.511176e+04 1.704993e-02 1.276e-02 3.232e+00
GE3 -2.511185e+04 -9.054298e-02 1.231e-03 3.298e+00
GE4 -2.511185e+04 -1.293560e-03 1.923e-03 3.277e+00
GE5 -2.511185e+04 -7.028176e-04 1.171e-03 3.294e+00
GE6 -2.511185e+04 2.954805e-04 9.963e-04 3.274e+00
GE7 -2.511185e+04 3.273638e-05 5.229e-04 3.229e+00
GE8 -2.511185e+04 -3.645170e-05 2.080e-04 3.244e+00
GE9 -2.511185e+04 1.642489e-05 2.044e-04 3.240e+00
GE10 -2.511185e+04 2.637428e-05 8.212e-05 3.342e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.852e-01 2.002e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 11
-------------------------------------------
DONE(9.557e+02 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511184e+04 0.000000e+00 1.022e-03 1.424e+01
GE2 -2.511184e+04 2.005819e-04 2.552e-03 3.258e+00
GE3 -2.511184e+04 -4.038308e-03 1.021e-03 3.258e+00
GE4 -2.511184e+04 -3.878337e-04 1.596e-04 3.344e+00
GE5 -2.511184e+04 -8.306884e-06 3.032e-04 3.298e+00
GE6 -2.511184e+04 7.114720e-05 9.693e-05 3.216e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.852e-01 2.001e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 12
-------------------------------------------
DONE(1.020e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511173e+04 0.000000e+00 5.402e-03 1.411e+01
GE2 -2.511171e+04 1.692736e-02 1.484e-02 3.287e+00
GE3 -2.511183e+04 -1.144752e-01 2.424e-03 3.353e+00
GE4 -2.511183e+04 -4.173679e-03 1.374e-03 3.231e+00
GE5 -2.511183e+04 -1.264840e-03 7.165e-04 3.253e+00
GE6 -2.511183e+04 -3.910133e-05 4.655e-04 3.214e+00
GE7 -2.511183e+04 4.997686e-05 4.037e-04 3.262e+00
GE8 -2.511183e+04 -1.149379e-06 2.200e-04 3.221e+00
GE9 -2.511183e+04 1.680431e-06 1.242e-04 3.284e+00
GE10 -2.511183e+04 4.469554e-06 5.130e-05 3.271e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.846e-01 1.997e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 13
-------------------------------------------
DONE(1.098e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511169e+04 0.000000e+00 5.681e-03 1.440e+01
GE2 -2.511169e+04 5.520196e-03 1.421e-02 3.324e+00
GE3 -2.511181e+04 -1.207179e-01 2.239e-03 3.329e+00
GE4 -2.511181e+04 -4.533073e-03 7.709e-04 3.311e+00
GE5 -2.511181e+04 -4.698135e-04 7.005e-04 3.269e+00
GE6 -2.511181e+04 -2.897758e-05 4.054e-04 3.300e+00
GE7 -2.511181e+04 5.393849e-05 3.201e-04 3.308e+00
GE8 -2.511181e+04 -7.190569e-06 1.554e-04 3.318e+00
GE9 -2.511181e+04 4.874389e-06 1.096e-04 3.342e+00
GE10 -2.511181e+04 4.319324e-06 4.203e-05 3.282e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.842e-01 1.995e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 14
-------------------------------------------
DONE(1.177e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511176e+04 0.000000e+00 2.519e-03 1.408e+01
GE2 -2.511175e+04 1.818197e-03 6.287e-03 3.244e+00
GE3 -2.511178e+04 -2.835488e-02 6.337e-04 3.218e+00
GE4 -2.511178e+04 -3.241519e-04 2.289e-04 3.250e+00
GE5 -2.511178e+04 -8.225882e-06 1.941e-04 3.239e+00
GE6 -2.511178e+04 1.325693e-05 1.862e-04 3.223e+00
GE7 -2.511178e+04 2.286369e-05 1.494e-04 3.217e+00
GE8 -2.511178e+04 1.716921e-05 5.752e-05 3.224e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.848e-01 1.998e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 15
-------------------------------------------
DONE(1.248e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511172e+04 0.000000e+00 2.722e-03 1.435e+01
GE2 -2.511171e+04 7.270926e-03 7.377e-03 3.260e+00
GE3 -2.511175e+04 -3.641411e-02 6.290e-04 3.232e+00
GE4 -2.511175e+04 -2.517177e-04 1.573e-04 3.236e+00
GE5 -2.511175e+04 -5.209650e-06 1.151e-04 3.253e+00
GE6 -2.511175e+04 7.782043e-06 1.063e-04 3.269e+00
GE7 -2.511175e+04 7.706228e-05 9.024e-05 3.348e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.853e-01 2.002e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 16
-------------------------------------------
DONE(1.316e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511163e+04 0.000000e+00 4.032e-03 1.428e+01
GE2 -2.511163e+04 -5.444345e-04 1.004e-02 3.251e+00
GE3 -2.511170e+04 -7.262161e-02 1.192e-03 3.251e+00
GE4 -2.511170e+04 -1.078686e-03 2.342e-04 3.347e+00
GE5 -2.511170e+04 -2.640850e-05 1.557e-04 3.232e+00
GE6 -2.511170e+04 1.513131e-07 1.077e-04 3.266e+00
GE7 -2.511170e+04 9.773332e-05 8.294e-05 3.262e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.849e-01 1.999e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 17
-------------------------------------------
DONE(1.383e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511164e+04 0.000000e+00 1.339e-03 1.434e+01
GE2 -2.511164e+04 1.863213e-04 3.031e-03 3.369e+00
GE3 -2.511165e+04 -7.322448e-03 3.050e-04 3.223e+00
GE4 -2.511165e+04 -1.580961e-05 9.223e-05 3.258e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.842e-01 1.995e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 18
-------------------------------------------
DONE(1.441e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511148e+04 0.000000e+00 4.606e-03 1.454e+01
GE2 -2.511150e+04 -1.285499e-02 9.536e-03 3.383e+00
GE3 -2.511159e+04 -9.148717e-02 7.370e-04 3.303e+00
GE4 -2.511159e+04 -4.969422e-04 1.561e-04 3.296e+00
GE5 -2.511159e+04 2.701909e-06 1.427e-04 3.318e+00
GE6 -2.511159e+04 1.295479e-06 1.106e-04 3.320e+00
GE7 -2.511159e+04 6.262190e-05 7.301e-05 3.262e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.844e-01 1.996e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 19
-------------------------------------------
DONE(1.509e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511147e+04 0.000000e+00 3.222e-03 1.444e+01
GE2 -2.511146e+04 4.390107e-03 8.230e-03 3.265e+00
GE3 -2.511151e+04 -5.172727e-02 6.026e-04 3.261e+00
GE4 -2.511151e+04 -4.358655e-04 1.426e-04 3.298e+00
GE5 -2.511151e+04 9.492484e-05 9.858e-05 3.288e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.850e-01 2.000e+03
------------------------------------------------------------------------------------------------
-------------------------------------------
STEP OF MOLECULAR DYNAMICS : 20
-------------------------------------------
DONE(1.570e+03 SEC) : INIT SCF
ITER ETOT(eV) EDIFF(eV) DRHO TIME(s)
GE1 -2.511135e+04 0.000000e+00 4.335e-03 1.445e+01
GE2 -2.511136e+04 -1.618388e-02 8.348e-03 3.293e+00
GE3 -2.511143e+04 -6.844721e-02 7.205e-04 3.326e+00
GE4 -2.511143e+04 -4.615375e-04 2.925e-04 3.254e+00
GE5 -2.511143e+04 1.542229e-05 2.092e-04 3.254e+00
GE6 -2.511143e+04 1.041771e-04 9.954e-05 3.240e+00
------------------------------------------------------------------------------------------------
Energy (Ry) Potential (Ry) Kinetic (Ry) Temperature (K)
-1.845e+03 -1.846e+03 2.850e-01 2.000e+03
------------------------------------------------------------------------------------------------
TIME STATISTICS
------------------------------------------------------------------------------------
CLASS_NAME NAME TIME(Sec) CALLS AVG(Sec) PER(%)
------------------------------------------------------------------------------------
total 1634.08 49 33.35 100.00
Driver reading 0.37 1 0.37 0.02
Input Init 0.36 1 0.36 0.02
Input_Conv Convert 0.00 1 0.00 0.00
Driver driver_line 1633.70 1 1633.70 99.98
UnitCell check_tau 0.00 1 0.00 0.00
PW_Basis_Sup setuptransform 0.02 1 0.02 0.00
PW_Basis_Sup distributeg 0.00 1 0.00 0.00
mymath heapsort 0.00 3 0.00 0.00
PW_Basis_K setuptransform 0.00 1 0.00 0.00
PW_Basis_K distributeg 0.00 1 0.00 0.00
PW_Basis setup_struc_factor 0.38 21 0.02 0.02
NOrbital_Lm extra_uniform 0.06 9 0.01 0.00
Mathzone_Add1 SplineD2 0.00 9 0.00 0.00
Mathzone_Add1 Cubic_Spline_Interpolation 0.00 9 0.00 0.00
Mathzone_Add1 Uni_Deriv_Phi 0.05 9 0.01 0.00
ppcell_vl init_vloc 0.01 1 0.01 0.00
Run_MD md_line 1633.13 1 1633.13 99.94
Nose_Hoover setup 102.89 1 102.89 6.30
MD_func force_virial 1632.22 21 77.72 99.89
ESolver_KS_LCAO Run 906.11 21 43.15 55.45
ESolver_KS_LCAO beforescf 7.77 21 0.37 0.48
ESolver_KS_LCAO beforesolver 1.99 21 0.09 0.12
ESolver_KS_LCAO set_matrix_grid 1.54 21 0.07 0.09
atom_arrange search 0.03 21 0.00 0.00
Grid_Technique init 1.50 21 0.07 0.09
Grid_BigCell grid_expansion_index 0.06 42 0.00 0.00
Record_adj for_2d 0.01 21 0.00 0.00
Grid_Driver Find_atom 0.04 4032 0.00 0.00
LCAO_Hamilt grid_prepare 0.00 21 0.00 0.00
OverlapNew initialize_SR 0.06 21 0.00 0.00
EkineticNew initialize_HR 0.02 21 0.00 0.00
NonlocalNew initialize_HR 0.21 21 0.01 0.01
Veff initialize_HR 0.02 21 0.00 0.00
LOC Alltoall 0.09 21 0.00 0.01
Charge set_rho_core 0.00 21 0.00 0.00
Charge atomic_rho 0.34 42 0.01 0.02
PW_Basis_Sup recip2real 2.22 1383 0.00 0.14
PW_Basis_Sup gathers_scatterp 0.86 1383 0.00 0.05
Potential init_pot 0.96 21 0.05 0.06
Potential update_from_charge 9.82 227 0.04 0.60
Potential cal_fixed_v 0.05 21 0.00 0.00
PotLocal cal_fixed_v 0.05 21 0.00 0.00
Potential cal_v_eff 9.74 227 0.04 0.60
H_Hartree_pw v_hartree 0.91 227 0.00 0.06
PW_Basis_Sup real2recip 2.26 1589 0.00 0.14
PW_Basis_Sup gatherp_scatters 0.82 1589 0.00 0.05
PotXC cal_v_eff 8.76 227 0.04 0.54
XC_Functional v_xc 8.72 227 0.04 0.53
Potential interpolate_vrs 0.03 227 0.00 0.00
H_Ewald_pw compute_ewald 0.05 21 0.00 0.00
Charge_Mixing init_mixing 0.00 21 0.00 0.00
HSolverLCAO solve 883.38 206 4.29 54.06
HamiltLCAO updateHk 541.51 206 2.63 33.14
OperatorLCAO init 394.87 618 0.64 24.16
OverlapNew calculate_SR 86.00 21 4.10 5.26
OverlapNew contributeHk 0.01 21 0.00 0.00
EkineticNew contributeHR 85.93 206 0.42 5.26
EkineticNew calculate_HR 85.92 21 4.09 5.26
NonlocalNew contributeHR 60.71 206 0.29 3.71
NonlocalNew calculate_HR 60.64 21 2.89 3.71
Veff contributeHR 308.67 206 1.50 18.89
Gint_interface cal_gint 753.77 433 1.74 46.13
Gint_interface cal_gint_vlocal 306.70 206 1.49 18.77
Gint_Tools cal_psir_ylm 127.53 741600 0.00 7.80
Gint_Gamma transfer_pvpR 1.90 206 0.01 0.12
OperatorLCAO contributeHk 0.11 206 0.00 0.01
HSolverLCAO hamiltSolvePsiK 14.49 206 0.07 0.89
OperatorLCAO get_hs_pointers 0.01 227 0.00 0.00
DiagoElpa elpa_solve 14.31 206 0.07 0.88
ElecStateLCAO psiToRho 327.38 206 1.59 20.03
ElecStateLCAO cal_dm_2d 2.38 206 0.01 0.15
elecstate cal_dm 3.12 227 0.01 0.19
psiMulPsiMpi pdgemm 3.03 227 0.01 0.19
DensityMatrix cal_DMR 0.15 206 0.00 0.01
Gint transfer_DMR 0.83 206 0.00 0.05
Gint_interface cal_gint_rho 324.01 206 1.57 19.83
Charge_Mixing get_drho 0.03 206 0.00 0.00
Charge mix_rho 1.13 185 0.01 0.07
Charge Pulay_mixing 0.12 185 0.00 0.01
ESolver_KS_LCAO out_deepks_labels 0.00 21 0.00 0.00
Force_Stress_LCAO getForceStress 726.10 21 34.58 44.44
Forces cal_force_loc 0.39 21 0.02 0.02
Forces cal_force_ew 0.34 21 0.02 0.02
Forces cal_force_cc 0.00 21 0.00 0.00
Forces cal_force_scc 0.48 21 0.02 0.03
Force_LCAO_gamma ftable_gamma 724.88 21 34.52 44.36
Force_LCAO_gamma allocate_gamma 419.68 21 19.98 25.68
LCAO_gen_fixedH b_NL_mu_new 173.74 21 8.27 10.63
Force_LCAO_gamma cal_foverlap 1.07 21 0.05 0.07
Force_LCAO_gamma cal_edm_2d 0.93 21 0.04 0.06
Force_LCAO_gamma cal_ftvnl_dphi 0.12 21 0.01 0.01
Force_LCAO_gamma cal_fvnl_dbeta_new 180.84 21 8.61 11.07
Force_LCAO_gamma cal_fvl_dphi 123.06 21 5.86 7.53
Gint_interface cal_gint_force 123.06 21 5.86 7.53
Gint_Tools cal_dpsir_ylm 69.81 37800 0.00 4.27
Nose_Hoover first_half 0.03 20 0.00 0.00
Nose_Hoover second_half 0.00 20 0.00 0.00
ModuleIO write_istate_info 0.01 1 0.01 0.00
------------------------------------------------------------------------------------
START Time : Mon Jun 3 08:12:42 2024
FINISH Time : Mon Jun 3 08:39:56 2024
TOTAL Time : 1634
SEE INFORMATION IN : OUT.Ti_nvt/
代码
文本
3.2 计算性能分析
根据输出日志，在2个CPU上运行，LAMMPS计算5000步的时间为786s，使用ABACUS第一性原理计算20步的时间为1634s。
使用深度势能平均每一步耗时0.157s，第一性原理每一步耗时81.7s，且系统大小只有前者的1/15。
代码
文本
4 参考文献
- Wang H et al. (2018). DeePMD-kit: A deep learning package for many-body potential energy representation and molecular dynamics. Comput. Phys. Commun. 228, 178-184.
- Wen T. et al. (2021). Specialising neural network potentials for accurate properties and application to the mechanical response of titanium. npj Computational Materials, 7, 1-11.
代码
文本
已赞1
推荐阅读
公开
以NaCl熔点为例学习DPMD梁大为
发布于 2024-06-01
公开
DPMD——分子动力学算法浅论邱翔宇
发布于 2024-06-02
1 赞