Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Local ISDF module (with k-sampling) for PBC mean-field calculation. #62

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
121 changes: 121 additions & 0 deletions examples/isdf/00-RHF_diamond.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,121 @@
import numpy as np
from pyscf import lib
from pyscf.gto.mole import *

from pyscf.isdf import isdf_tools_cell
from pyscf.isdf import isdf_local_k
from pyscf.isdf import isdf_jk
from pyscf.isdf import isdf_local

from pyscf.lib.parameters import BOHR

MOL_STRUCTURE = '''
C 0. 0. 0.
C 0.8917 0.8917 0.8917
C 1.7834 1.7834 0.
C 2.6751 2.6751 0.8917
C 1.7834 0. 1.7834
C 2.6751 0.8917 2.6751
C 0. 1.7834 1.7834
C 0.8917 2.6751 2.6751
'''

#### NOTE: a full tests on combinations of parameters ####

C_ARRAY = [15, 15, 20, 25, 30, 30]
RELA_CUTOFF = [3e-2, 1e-2, 3e-3, 1e-3, 3e-4, 1e-4]
SuperCell_ARRAY = [
# [1, 1, 1],
[1, 1, 2],
# [1, 2, 2],
# [2, 2, 2],
# [3, 3, 3],
# [2, 4, 4],
# [3, 4, 4],
# [5, 5, 5],
# [6, 6, 6],
# [1, 1, 4],
# [1, 1, 8],
# [1, 1, 16],
# [1, 1, 32],
]


Ke_CUTOFF = [70]
boxlen = 3.5668
Basis = ['gth-dzvp']

PARTITION = [
[[0,1],[2,3],[4,5],[6,7]],
[[0,1,2,3],[4,5,6,7]],
[[0,1,2,3,4,5,6,7]],
[[0],[1],[2],[3],[4],[5],[6],[7]],
]

if __name__ == '__main__':

boxlen = 3.57371000
prim_a = np.array([[boxlen,0.0,0.0],[0.0,boxlen,0.0],[0.0,0.0,boxlen]])
atm = [
['C', (0. , 0. , 0. )],
['C', (0.8934275 , 0.8934275 , 0.8934275)],
['C', (1.786855 , 1.786855 , 0. )],
['C', (2.6802825 , 2.6802825 , 0.8934275)],
['C', (1.786855 , 0. , 1.786855)],
['C', (2.6802825 , 0.8934275 , 2.6802825)],
['C', (0. , 1.786855 , 1.786855)],
['C', (0.8934275 , 2.6802825 , 2.6802825)],
]

for supercell in SuperCell_ARRAY:
ke_cutoff = Ke_CUTOFF[0]
for partition in PARTITION: ## test different partition of atoms
for basis in Basis:
for c, rela_cutoff in zip(C_ARRAY, RELA_CUTOFF):
# for c in C_ARRAY:
print('--------------------------------------------')
print('C = %.2e, supercell = %s, kc_cutoff = %d, basis = %s, partition = %s' % (
c, str(supercell), ke_cutoff, basis, partition))

prim_cell = isdf_tools_cell.build_supercell(atm, prim_a, Ls = [1,1,1], ke_cutoff=ke_cutoff, basis=basis, pseudo="gth-pade", verbose=4)
prim_mesh = prim_cell.mesh
print("prim_mesh = ", prim_mesh)

mesh = [supercell[0] * prim_mesh[0], supercell[1] * prim_mesh[1], supercell[2] * prim_mesh[2]]
mesh = np.array(mesh, dtype=np.int32)

cell, supercell_group = isdf_tools_cell.build_supercell_with_partition(atm, prim_a, partition=partition, Ls = supercell, ke_cutoff=ke_cutoff, mesh=mesh, basis=basis, pseudo="gth-pade", verbose=4)

cell.incore_anyway = False
cell.max_memory = 200 # force to call with_df.get_jk

t1 = (lib.logger.process_clock(),lib.logger.perf_counter())

pbc_isdf_info = isdf_local.PBC_ISDF_Info_Quad(cell, with_robust_fitting=True, direct=False, rela_cutoff_QRCP=rela_cutoff)
pbc_isdf_info.build_IP_local(c=c, group=supercell_group, Ls=[supercell[0]*4, supercell[1]*4, supercell[2]*4])
print("pbc_isdf_info.naux = ", pbc_isdf_info.naux)
print("effective c = ", float(pbc_isdf_info.naux) / pbc_isdf_info.nao)
pbc_isdf_info.build_auxiliary_Coulomb()

t2 = (lib.logger.process_clock(), lib.logger.perf_counter())

print(isdf_jk._benchmark_time(t1, t2, 'build_isdf', pbc_isdf_info))

# for bunch_size in BUNCHSIZE_IO:
### perform scf ###

from pyscf.pbc import scf

t1 = (lib.logger.process_clock(), lib.logger.perf_counter())
mf = scf.RHF(cell)
mf.with_df = pbc_isdf_info
mf.max_cycle = 32
mf.conv_tol = 1e-7
pbc_isdf_info.direct_scf = mf.direct_scf
mf.kernel()
t2 = (lib.logger.process_clock(), lib.logger.perf_counter())
print(isdf_jk._benchmark_time(t1, t2, 'scf_isdf', pbc_isdf_info))

del mf
del pbc_isdf_info
exit(1)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there a reason why you exit(1) in all of the examples? This signifies an error has occurred which should not be the case if it reaches the end of the code.

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for pointing out this issue. It is added just because I do not want run all the for-loops when I am testing the code. It should not appear here and I will comment this out.

115 changes: 115 additions & 0 deletions examples/isdf/01-KRHF_TiO2.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,115 @@
from functools import reduce
import numpy as np
from pyscf import lib
import pyscf.pbc.gto as pbcgto
from pyscf.pbc.gto import Cell
from pyscf.pbc import tools
from pyscf.pbc.lib.kpts import KPoints
from pyscf.pbc.lib.kpts_helper import is_zero, gamma_point, member
from pyscf.gto.mole import *

from pyscf.isdf import isdf_tools_cell
from pyscf.isdf import isdf_local_k
from pyscf.isdf import isdf_jk

MOL_STRUCTURE = '''
Ti 2.3246330643 2.3246330643 1.4853414945
Ti 0.0000000000 0.0000000000 -0.0000000000
O 0.9065353261 3.7427308025 1.4853414945
O 3.7427308025 0.9065353261 1.4853414945
O 1.4180977382 1.4180977382 0.0000000000
O 3.2311683903 3.2311683903 0.0000000000
'''

atm = [
['Ti',(2.3246330643,2.3246330643, 1.4853414945)],
['Ti',(0.0000000000,0.0000000000, 0.0000000000)],
['O ',(0.9065353261,3.7427308025, 1.4853414945)],
['O ',(3.7427308025,0.9065353261, 1.4853414945)],
['O ',(1.4180977382,1.4180977382, 0.0000000000)],
['O ',(3.2311683903,3.2311683903, 0.0000000000)],
]
boxlen = [4.6492659759,4.6492659759,2.9706828877]

C_ARRAY = [15,20,25,30] ## if rela_cutoff_QRCP is set, then c is used to when performing random projection, which can be relative large.
RELA_QR = [1e-2,1e-3,2e-4,1e-4]
SuperCell_ARRAY = [
# [1,1,1],
[2,2,2],
[3,3,3],
[4,4,4],
[5,5,5],
[6,6,6],
]
Ke_CUTOFF = [128, 192]

Basis = ['gth-cc-tzvp-Ye']

prim_partition = [[0],[1],[2],[3],[4],[5]]

if __name__ == '__main__':

prim_a = np.array([[boxlen[0],0.0,0.0],[0.0,boxlen[1],0.0],[0.0,0.0,boxlen[2]]])
pseudo = 'gth-hf-rev'

for supercell in SuperCell_ARRAY:
for basis in Basis:
for ke_cutoff in Ke_CUTOFF:

DM_CACHED = None

from pyscf.gto.basis import parse_nwchem
fbas="basis.dat"
atms = ['O', 'Ti']
basis = {atm:parse_nwchem.load(fbas, atm) for atm in atms}
print("basis = ", basis)


prim_cell = isdf_tools_cell.build_supercell(atm, prim_a, Ls = [1,1,1], ke_cutoff=ke_cutoff, basis=basis, pseudo=pseudo, spin=0, verbose=10)
cell = prim_cell

### perform scf ###

from pyscf.pbc import scf, dft
from pyscf.pbc.dft import multigrid

nk = supercell
kpts = cell.make_kpts(nk)

for c,rela_qr in list(zip(C_ARRAY,RELA_QR)):

print('--------------------------------------------')
print('C = %d, QR=%f, supercell = %s, kc_cutoff = %d, basis = %s' % (c, rela_qr, str(supercell), ke_cutoff, basis))

### create the isdf object ###

t1 = (lib.logger.process_clock(), lib.logger.perf_counter())
pbc_isdf_info = isdf_local_k.PBC_ISDF_Info_Quad_K(cell,
kmesh=nk,
with_robust_fitting=True,
rela_cutoff_QRCP=rela_qr,
direct=True,
limited_memory=True,
build_K_bunchsize=128, ## NOTE:control the memory cost in building K
# use_occ_RI_K=False
)
pbc_isdf_info.verbose = 10
pbc_isdf_info.build_IP_local(c=c, m=5, group=prim_partition)
print("effective c = ", float(pbc_isdf_info.naux) / pbc_isdf_info.nao)
t2 = (lib.logger.process_clock(), lib.logger.perf_counter())
print(isdf_jk._benchmark_time(t1, t2, 'build ISDF', pbc_isdf_info))

t1 = (lib.logger.process_clock(), lib.logger.perf_counter())
mf = scf.KRHF(cell, kpts)
mf.with_df = pbc_isdf_info
mf.max_cycle = 100
mf.conv_tol = 1e-8
mf.conv_tol_grad = 1e-3
if DM_CACHED is not None:
mf.kernel(DM_CACHED)
else:
mf.kernel()
t2 = (lib.logger.process_clock(), lib.logger.perf_counter())

print(isdf_jk._benchmark_time(t1, t2, 'RHF_bench', mf))
DM_CACHED = mf.make_rdm1()
139 changes: 139 additions & 0 deletions examples/isdf/02-UHF_CCO.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,139 @@
import numpy as np
from pyscf import lib
from pyscf.gto.mole import *

from pyscf.isdf import isdf_tools_cell
from pyscf.isdf import isdf_local_k
from pyscf.isdf import isdf_jk
from pyscf.isdf import isdf_local

from pyscf.lib.parameters import BOHR

#### NOTE: a full tests on combinations of parameters ####

prim_a = np.array(
[[14.572056092, 0.000000000, 0.000000000],
[0.000000000, 14.572056092, 0.000000000],
[0.000000000, 0.000000000, 6.010273939],]) * BOHR
atm = [
['Cu', (1.927800, 1.927800, 1.590250)],
['Cu', (5.783400, 5.783400, 1.590250)],
['Cu', (1.927800, 5.783400, 1.590250)],
['Cu', (5.783400, 1.927800, 1.590250)],
['O', (1.927800, 3.855600, 1.590250)],
['O', (3.855600, 5.783400, 1.590250)],
['O', (5.783400, 3.855600, 1.590250)],
['O', (3.855600, 1.927800, 1.590250)],
['O', (0.000000, 1.927800, 1.590250)],
['O', (1.927800, 7.711200, 1.590250)],
['O', (7.711200, 5.783400, 1.590250)],
['O', (5.783400, 0.000000, 1.590250)],
['Ca', (0.000000, 0.000000, 0.000000)],
['Ca', (3.855600, 3.855600, 0.000000)],
['Ca', (7.711200, 3.855600, 0.000000)],
['Ca', (3.855600, 7.711200, 0.000000)],
]

C_ARRAY = [25, 30, 35]
RELA_CUTOFF = [1e-3, 3e-4, 1e-4]
SuperCell_ARRAY = [
[1, 1, 1],
]
Ke_CUTOFF = [256]
Basis = ['gth-dzvp']

PARTITION = [
[[0], [1], [2], [3],
[4], [5], [6], [7],
[8], [9], [10], [11],
[12], [13], [14], [15]]
]

if __name__ == '__main__':

for supercell in SuperCell_ARRAY:
ke_cutoff = Ke_CUTOFF[0]
for partition in PARTITION: ## test different partition of atoms
for _basis_ in Basis:

DM_CACHED = None

from pyscf.gto.basis import parse_nwchem
fbas="basis2.dat"
atms = ['O', 'Cu', "Ca"]
basis = {atm:parse_nwchem.load(fbas, atm) for atm in atms}
# print("basis = ", basis)

pseudo = {'Cu': 'gth-pbe-q19', 'O': 'gth-pbe', 'Ca': 'gth-pbe'}

prim_cell = isdf_tools_cell.build_supercell(atm, prim_a, Ls = [1,1,1], ke_cutoff=ke_cutoff, basis=basis, pseudo=pseudo, verbose=4)
prim_mesh = prim_cell.mesh
# print("prim_mesh = ", prim_mesh)

mesh = [supercell[0] * prim_mesh[0], supercell[1] * prim_mesh[1], supercell[2] * prim_mesh[2]]
mesh = np.array(mesh, dtype=np.int32)

cell, supercell_group = isdf_tools_cell.build_supercell_with_partition(atm, prim_a,
partition = partition,
Ls = supercell,
ke_cutoff = ke_cutoff,
mesh = mesh,
basis = basis,
pseudo = pseudo,
verbose = 4)

cell.incore_anyway = False
cell.max_memory = 200 # force to call with_df.get_jk

for c, rela_cutoff in zip(C_ARRAY, RELA_CUTOFF):

print('--------------------------------------------')
print('C = %.2e, supercell = %s, kc_cutoff = %d, basis = %s, partition = %s' % (
c, str(supercell), ke_cutoff, basis, partition))

t1 = (lib.logger.process_clock(),lib.logger.perf_counter())

pbc_isdf_info = isdf_local.PBC_ISDF_Info_Quad(cell, with_robust_fitting=True,
direct=True,
rela_cutoff_QRCP=rela_cutoff,
limited_memory=True, build_K_bunchsize=56
)
pbc_isdf_info.build_IP_local(c=c, group=supercell_group)
print("pbc_isdf_info.naux = ", pbc_isdf_info.naux)
print("effective c = ", float(pbc_isdf_info.naux) / pbc_isdf_info.nao)
pbc_isdf_info.build_auxiliary_Coulomb()

t2 = (lib.logger.process_clock(), lib.logger.perf_counter())

print(isdf_jk._benchmark_time(t1, t2, 'build_isdf', pbc_isdf_info))

# for bunch_size in BUNCHSIZE_IO:
### perform scf ###

from pyscf.pbc import scf

t1 = (lib.logger.process_clock(), lib.logger.perf_counter())
mf = scf.UHF(cell)
mf.with_df = pbc_isdf_info
mf.max_cycle = 64
mf.conv_tol = 1e-7
pbc_isdf_info.direct_scf = mf.direct_scf
if DM_CACHED is not None:
mf.kernel(DM_CACHED)
else:
mf.kernel()
t2 = (lib.logger.process_clock(), lib.logger.perf_counter())
print(isdf_jk._benchmark_time(t1, t2, 'scf_isdf', pbc_isdf_info))

del mf
del pbc_isdf_info

### GDF benchmark ###

mf = scf.UHF(cell).density_fit()
mf.max_cycle = 64
mf.conv_tol = 1e-7
# pbc_isdf_info.direct_scf = mf.direct_scf
mf.kernel(DM_CACHED)

exit(1)
Loading
Loading