aiida-pyscf

一个 AiiDA 插件，用于基于Python的化学模拟框架（PySCF）。

1. 安装
2. 需求
3. 设置
4. 示例
   • 平均场计算
   • 自定义结构
   • 优化几何形状
   • 将哈密顿量写入FCIDUMP文件
   • 将轨道写入CUBE文件
   • 重启未收敛的计算
   • 自动错误恢复
   • 序列化的模型

安装

推荐的安装方法是通过 pip

pip install aiida-pyscf

需求

要使用 aiida-pyscf，需要一个配置好的AiiDA配置文件。请参阅 aiida-core的文档 以获取详细说明。

设置

要通过AiiDA使用 aiida-pyscf 插件运行PySCF计算，需要在计算机上配置PySCF的运行位置。请参阅 aiida-core的文档 以获取详细说明。

然后需要配置PySCF代码。以下YAML配置文件可以作为起点

label: pyscf
description: PySCF
computer: localhost
filepath_executable: python
default_calc_job_plugin: pyscf.base
use_double_quotes: false
with_mpi: false
prepend_text: ''
append_text: ''

将内容写入名为 pyscf.yml 的文件中，确保将 computer 的值更新为上一步配置的计算机标签。要配置代码，请执行

verdi code create core.code.installed --config pyscf.yml -n

现在应该已经创建了一个带有标签 pyscf 的代码，该代码将在下面的示例中使用。

示例

平均场计算

默认的计算是执行平均场计算。至少，应该定义结构和平均场方法。

from ase.build import molecule
from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code

builder = load_code('pyscf').get_builder()
builder.structure = StructureData(ase=molecule('H2O'))
builder.parameters = Dict({'mean_field': {'method': 'RHF'}})
results, node = run.get_node(builder)

该代码在水分子的几何结构上运行 Hartree-Fock 计算。

主要结果存储在 parameters 输出中，默认包含计算得到的 total_energy 和 forces，分子轨道的详细信息以及一些计时信息。

print(results['parameters'].get_dict())
{
    'mean_field': {
        'forces': [
            [-6.4898366104394e-16, 3.0329042995656e-15, 2.2269765466236],
            [1.122487932593e-14, 0.64803103141326, -1.1134882733107],
            [-1.0575895664886e-14, -0.64803103141331, -1.1134882733108]
        ],
        'forces_units': 'eV/Å',
        'molecular_orbitals': {
            'labels': [
                '0 O 1s',
                '0 O 2s',
                '0 O 2px',
                '0 O 2py',
                '0 O 2pz',
                '1 H 1s',
                '2 H 1s'
            ],
            'energies': [
                -550.86280025028,
                -34.375426862456,
                -16.629598134599,
                -12.323304634736,
                -10.637428057751,
                16.200273277782,
                19.796075801491
            ],
            'occupations': [2.0, 2.0, 2.0, 2.0, 2.0, 0.0, 0.0]
        },
        'total_energy': -2039.8853743664,
        'total_energy_units': 'eV',
    },
    'timings': {
        'total': 1.3238215579768, 'mean_field': 0.47364449803717
    },
}

自定义结构

结构的几何形状通过 structure 输入完全定义，该输入由 StructureData 节点提供。任何其他属性，例如电荷和要使用的基组，都可以通过 parameters 输入中的 structure 字典进行指定。

from ase.build import molecule
from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code

builder = load_code('pyscf').get_builder()
builder.structure = StructureData(ase=molecule('H2O'))
builder.parameters = Dict({
    'mean_field': {'method': 'RHF'},
    'structure': {
        'basis ': 'sto-3g',
        'charge': 0,
    }
})
results, node = run.get_node(builder)

用于定义结构的 pyscf.gto.Mole 类 的任何属性都可以通过 structure 字典设置，除了由插件根据 StructureData 输入自动设置的 atom 和 unit 属性。

优化几何形状

可以通过在输入 parameters 中指定 optimizer 字典来优化几何形状。必须指定 solver，目前支持的求解器有 geometric 和 berny。convergence_parameters 接受所选求解器的参数（有关详细信息，请参阅 PySCF 文档）。

from ase.build import molecule
from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code

builder = load_code('pyscf').get_builder()
builder.structure = StructureData(ase=molecule('H2O'))
builder.parameters = Dict({
    'mean_field': {'method': 'RHF'},
    'optimizer': {
        'solver': 'geometric',
        'convergence_parameters': {
            'convergence_energy': 1e-6,  # Eh
            'convergence_grms': 3e-4,    # Eh/Bohr
            'convergence_gmax': 4.5e-4,  # Eh/Bohr
            'convergence_drms': 1.2e-3,  # Angstrom
            'convergence_dmax': 1.8e-3,  # Angstrom
        }
    }
})
results, node = run.get_node(builder)

优化后的结构以 StructureData 的形式返回，标签为 structure。几何优化轨迹中每个框架的结构和能量以 TrajectoryData 的形式存储，标签为 trajectory。总能量可以通过以下方式检索：

results['trajectory'].get_array('energies')

定位轨道

要计算定位轨道，请在 parameters.localize_orbitals.method 输入中指定所需的方法。

from ase.build import molecule
from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code

builder = load_code('pyscf').get_builder()
builder.structure = StructureData(ase=molecule('H2O'))
builder.parameters = Dict({
    'mean_field': {'method': 'RHF'},
    'localize_orbitals': {'method': 'ibo'}
})
results, node = run.get_node(builder)

支持以下方法：boys、cholesky、edmiston、iao、ibo、lowdin、nao、orth、pipek、vvo。有关更多信息，请参阅 PySCF 文档。

计算 Hessian

为了计算 Hessian，请在 parameters 输入中将 hessian 键指定为空字典。

from ase.build import molecule
from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code

builder = load_code('pyscf').get_builder()
builder.structure = StructureData(ase=molecule('H2O'))
builder.parameters = Dict({
    'mean_field': {'method': 'RHF'},
    'hessian': {}
})
results, node = run.get_node(builder)

计算得到的 Hessian 将作为带有 hessian 链接标签的 ArrayData 节点附加。使用 node.outputs.hessian.get_array('hessian') 将计算得到的 Hessian 作为 numpy 数组检索以进行进一步处理。

将哈密顿量写入FCIDUMP文件

要指示计算将 Hamiltonian 的表示写入 FCIDUMP 文件，请将 fcidump 字典添加到 parameters 输入中。

from ase.build import molecule
from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code

builder = load_code('pyscf').get_builder()
builder.structure = StructureData(ase=molecule('N2'))
builder.parameters = Dict({
    'mean_field': {'method': 'RHF'},
    'fcidump': {
        'active_spaces': [[5, 6, 8, 9]],
        'occupations': [[1, 1, 1, 1]]
    }
})
results, node = run.get_node(builder)

active_spaces 和 occupations 键是必需的，并且每个都接受一个整数列表的列表。列表中的每个元素都会为相应的活跃空间和轨道的占据生成一个 FCIDUMP 文件。active_spaces 和 occupations 数组的形状必须相同。

生成的 FCIDUMP 文件作为 fcidump 命名空间中的 SinglefileData 输出节点附加，标签由列表中相应活跃空间的索引确定。

print(results['fcidump']['active_space_0'].get_content())
 &FCI NORB=   4,NELEC= 4,MS2=0,
  ORBSYM=1,1,1,1,
  ISYM=1,
 &END
  0.5832127121682998       1    1    1    1
  0.5359642500498074       1    1    2    2
 -2.942091015256668e-15    1    1    3    2
  0.5381290185905914       1    1    3    3
 -3.782672959584676e-15    1    1    4    1
  ...

生成 CUBE 文件

pyscf.tools.cubegen 模块提供用于计算系统各种属性并将其写入 CUBE 文件的函数。目前，PyscfCalculation 插件支持以下计算：

• 分子轨道
• 电荷密度
• 分子静电势

要指示计算将任何这些量的表示写入 CUBE 文件，请将 cubegen 字典添加到 parameters 输入中。

from ase.build import molecule
from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code

builder = load_code('pyscf').get_builder()
builder.structure = StructureData(ase=molecule('N2'))
builder.parameters = Dict({
    'mean_field': {'method': 'RHF'},
    'cubegen': {
        'orbitals: {
            'indices': [5, 6],
            'parameters': {
                'nx': 40,
                'ny': 40,
                'nz': 40,
            }
        },
        'density': {
            'parameters': {
                'resolution': 300,
            }
        },
        'mep': {
            'parameters': {
                'resolution': 300,
            }
        }
    }
})
results, node = run.get_node(builder)

在 orbitals 子字典中必须指定 indices 键，它接受一个整数列表，表示应写入文件的分子轨道的索引。可以在 parameters 子字典中提供其他参数（有关详细信息，请参阅PySCF 文档）。对于 density 和 mep 字典，parameters 子字典是可选的。要计算电荷密度和分子静电势，分别需要为 density 和 mep 键提供一个空字典。

生成的 CUBE 文件作为 SinglefileData 输出节点附加在 cubegen 命名空间中，具有 orbitals、density 和 mep 子命名空间。对于 orbitals 子命名空间，标签由相应的分子轨道索引确定

print(results['cubegen']['orbitals']['mo_5'].get_content())
Orbital value in real space (1/Bohr^3)
PySCF Version: 2.1.1  Date: Sun Apr  2 15:59:19 2023
    2   -3.000000   -3.000000   -4.067676
   40    0.153846    0.000000    0.000000
   40    0.000000    0.153846    0.000000
   40    0.000000    0.000000    0.208599
    7    0.000000    0.000000    0.000000    1.067676
    7    0.000000    0.000000    0.000000   -1.067676
 -1.10860E-04 -1.56874E-04 -2.16660E-04 -2.92099E-04 -3.84499E-04 -4.94299E-04
 -6.20809E-04 -7.62048E-04 -9.14724E-04 -1.07439E-03 -1.23579E-03 -1.39331E-03
  ...

警告 已知 PySCF 在使用 DHF、DKS、GHF 和 GKS 引用时计算 MEP 会失败。

重启未收敛的计算

插件将自动指导 PySCF 编写检查点文件。如果计算没有收敛，它将以退出状态 410 结束，并将检查点文件作为 SinglefileData 附加到 checkpoint 输出节点。然后可以将此节点作为输入传递给新的计算，从检查点重新开始

failed_calculation = load_node(IDENTIFIER)
builder = failed_calculation.get_builder_restart()
builder.checkpoint = failed_calculation.outputs.checkpoint
submit(builder)

插件将把失败的计算的检查点文件写入工作目录，以便 PySCF 可以从这里开始。

后处理

PyscfCalculation 插件不支持所有 PySCF 功能；它旨在支持大多数计算密集型的功能，如在此例中，能够将计算作为远程计算资源上的 calcjob 脱载。大多数后处理实用程序计算成本低廉，并且由于 API 是 Python，它们可以直接在 AiiDA 工作流程中以 calcfunction 的形式调用。许多 PySCF 实用程序需要系统的 model 作为参数，其中 model 是 PySCF 中使用的主体对象，即以下代码中分配给 mean_field 变量的对象

from pyscf import scf
mean_field = scf.RHF(..)
mean_field.kernel()

kernel 方法通常计算成本很高，但它的结果（存储在模型对象中）在 PyscfCalculation 结束时丢失，因为 calcjob 的 Python 解释器关闭，并且 mean_field 对象不再存在。这会迫使后处理代码从头开始重建模型并重新运行昂贵的内核。因此，PyscfCalculation 将计算出的 PySCF 模型序列化，并将其存储为带有链接标签 model 的 PickledData 输出节点，在归因图中。这允许在另一个 Python 解释器中重新创建模型，并使其可用于后处理

from pyscf.hessian import thermo
node = load_node()  # Load the completed `PyscfCalculation`
mean_field = node.outputs.model.load()  # Reconstruct the model by calling the `load()` method
hessian = mean_field.Hessian().kernel()
freq_info = thermo.harmonic_analysis(mean_field.mol, hessian)

自动错误恢复

PySCF 计算可能无法以预期结果结束的原因有很多。例如，自洽场循环无法收敛或作业因调度程序杀死而超出了请求的 walltime。PyscfBaseWorkChain 被设计成在可能的情况下自动从这些错误中恢复。它是围绕 PyscfCalculation 插件的一个简单包装，如果前一次迭代失败，则自动重新启动新的 PyscfCalculation。启动 PyscfBaseWorkChain 几乎与直接启动 PyscfCalculation 一样；输入只需在 pyscf 命名空间内“嵌套”即可

from aiida.engine import run
from aiida.orm import Dict, StructureData, load_code, load_node
from aiida_pyscf.workflows.base import PyscfBaseWorkChain
from ase.build import molecule

builder = PyscfBaseWorkChain.get_builder()
builder.pyscf.code = load_code('pyscf')
builder.pyscf.structure = StructureData(ase=molecule('H2O'))
builder.pyscf.parameters = Dict({
    'mean_field': {
        'method': 'RHF',
        'max_cycle': 3,
    }
})
results, node = run.get_node(builder)

在此示例中，我们故意将自洽场循环的最大迭代次数设置为 3（《'mean_field.max_cycle' = 3》），这将导致第一次迭代无法达到收敛。PyscfBaseWorkChain 检测到错误，由 PyscfCalculation 上的退出状态 410 标识，并自动从保存的检查点重新启动计算。经过三次迭代，计算收敛

$ verdi process status IDENTIFIER
PyscfBaseWorkChain<30126> Finished [0] [2:results]
    ├── PyscfCalculation<30127> Finished [410]
    ├── PyscfCalculation<30132> Finished [410]
    └── PyscfCalculation<30137> Finished [0]

以下错误模式目前由 PyscfBaseWorkChain 处理：

• 120：超出运行时间：如果可用，计算将从最后一个检查点重新启动，否则工作链将终止
• 140：节点故障：计算将从最后一个检查点重新启动
• 410：电子收敛未实现：计算将从最后一个检查点重新启动
• 500：离子收敛未实现：几何优化未收敛，计算将从最后一个检查点和结构重新启动

序列化的模型

PyscfCalculation 的主要目标是解决给定结构的平均场问题。这一步骤的结果，通常计算成本较高，存储在主脚本中的 mean_field_run 变量中

mean_field = scf.RHF(structure)
density_matrix = mean_field.from_chk('restart.chk')
mean_field_run = mean_field.run(density_matrix)

可以使用 mean_field_run 对象在 PySCF 中实现许多后续后处理操作。为了保持 PyscfCalculation 接口简单，并非所有这些功能都得到支持。然而，一旦计算作业完成，mean_field_run 变量就会丢失，无法再用于进一步处理。

作为解决方案，PyscfCalculation 将 mean_field_run 对象 "pickle" 并将其附加到计算的 model 输出。可以“unpickle” model 输出节点以恢复原始的 mean_field_run 对象，以便用于进一步处理

from aiida.engine import run
inputs = {}
results, node = run.get_node(PyscfCalculation, **inputs)
mean_field = node.outputs.model.load()
print(mean_field.e_tot)

警告 对于某些情况，计算可能无法 pickle 模型并引发异常。在这种情况下，可以将 PyscfCalculation 的 pickle_model 输入设置为 False。

贡献

此项目欢迎贡献和建议。大多数贡献需要您同意贡献者许可协议（CLA），声明您有权，并且实际上确实授予我们使用您的贡献的权利。有关详细信息，请访问 https://cla.opensource.microsoft.com。

当您提交拉取请求时，CLA 机器人将自动确定您是否需要提供 CLA，并相应地装饰 PR（例如，状态检查，评论）。只需遵循机器人提供的说明即可。您只需在整个使用我们的 CLA 的所有存储库中执行此操作一次。

此项目已采用 Microsoft 开源行为准则。有关更多信息，请参阅 行为准则常见问题解答 或联系 opencode@microsoft.com 以提出任何额外的问题或评论。

商标

此项目可能包含项目、产品或服务的商标或徽标。Microsoft 商标或徽标的授权使用必须遵守并遵循 Microsoft 的商标和品牌指南。在此项目的修改版本中使用 Microsoft 商标或徽标不得引起混淆或暗示 Microsoft 赞助。第三方商标或徽标的任何使用均受这些第三方政策的约束。