diff --git a/gpu4pyscf/df/int3c2e.py b/gpu4pyscf/df/int3c2e.py
index 523dc5e4..d4171d25 100644
--- a/gpu4pyscf/df/int3c2e.py
+++ b/gpu4pyscf/df/int3c2e.py
@@ -1453,118 +1453,6 @@ def get_int3c2e(mol, auxmol=None, auxbasis='weigend+etb', direct_scf_tol=1e-13,
return int3c.transpose([2,1,0])
-def get_int2c2e_sorted(mol, auxmol, intopt=None, direct_scf_tol=1e-13, aosym=None, omega=None, stream=None):
- '''
- Generated int2c2e consistent with pyscf
- '''
- if omega is None: omega = 0.0
- if stream is None: stream = cupy.cuda.get_current_stream()
- if intopt is None:
- intopt = VHFOpt(mol, auxmol, 'int2e')
- intopt.build(direct_scf_tol, diag_block_with_triu=True, aosym=False)
- naux = auxmol.nao
- rows, cols = np.tril_indices(naux)
-
- nbins = 1
-
- nao_cart = intopt.mol.nao
- naux_cart = intopt.auxmol.nao
- norb_cart = nao_cart + naux_cart + 1
-
- int2c = cupy.zeros([naux_cart, naux_cart], order='F')
- ao_offsets = np.array([nao_cart+1, nao_cart, nao_cart+1, nao_cart], dtype=np.int32)
- strides = np.array([1, naux_cart, naux_cart, naux_cart*naux_cart], dtype=np.int32)
- for k_id, log_q_k in enumerate(intopt.aux_log_qs):
- bins_locs_k = _make_s_index_offsets(log_q_k, nbins)
- cp_k_id = k_id + len(intopt.log_qs)
- for l_id, log_q_l in enumerate(intopt.aux_log_qs):
- if k_id > l_id: continue
- bins_locs_l = _make_s_index_offsets(log_q_l, nbins)
- cp_l_id = l_id + len(intopt.log_qs)
- err = libgint.GINTfill_int2e(
- ctypes.cast(stream.ptr, ctypes.c_void_p),
- intopt.bpcache,
- ctypes.cast(int2c.data.ptr, ctypes.c_void_p),
- ctypes.c_int(norb_cart),
- strides.ctypes.data_as(ctypes.c_void_p),
- ao_offsets.ctypes.data_as(ctypes.c_void_p),
- bins_locs_k.ctypes.data_as(ctypes.c_void_p),
- bins_locs_l.ctypes.data_as(ctypes.c_void_p),
- ctypes.c_int(nbins),
- ctypes.c_int(cp_k_id),
- ctypes.c_int(cp_l_id),
- ctypes.c_double(omega))
-
- if err != 0:
- raise RuntimeError("int2c2e failed\n")
-
- int2c[rows, cols] = int2c[cols, rows]
- if not mol.cart:
- coeff = intopt.aux_cart2sph
- int2c = coeff.T @ int2c @ coeff
-
- return int2c
-
-def get_int2c2e_ip_sorted(mol, auxmol, intopt=None, direct_scf_tol=1e-13, intor=None, aosym=None, stream=None):
- '''
- TODO: WIP
- '''
- if stream is None: stream = cupy.cuda.get_current_stream()
- if intopt is None:
- intopt = VHFOpt(mol, auxmol, 'int2e')
- intopt.build(direct_scf_tol, diag_block_with_triu=True, aosym=False)
-
- nbins = 1
-
- nao_cart = intopt.mol.nao
- naux_cart = intopt.auxmol.nao
- norb_cart = nao_cart + naux_cart + 1
- rows, cols = np.tril_indices(naux_cart)
-
- int2c = cupy.zeros([naux_cart, naux_cart], order='F')
- ao_offsets = np.array([nao_cart+1, nao_cart, nao_cart+1, nao_cart], dtype=np.int32)
- strides = np.array([1, naux_cart, naux_cart, naux_cart*naux_cart], dtype=np.int32)
- for k_id, log_q_k in enumerate(intopt.aux_log_qs):
- bins_locs_k = _make_s_index_offsets(log_q_k, nbins)
- cp_k_id = k_id + len(intopt.log_qs)
- for l_id, log_q_l in enumerate(intopt.aux_log_qs):
- if k_id > l_id: continue
- bins_locs_l = _make_s_index_offsets(log_q_l, nbins)
- cp_l_id = l_id + len(intopt.log_qs)
- err = libgint.GINTfill_int2e(
- ctypes.cast(stream.ptr, ctypes.c_void_p),
- intopt.bpcache,
- ctypes.cast(int2c.data.ptr, ctypes.c_void_p),
- ctypes.c_int(norb_cart),
- strides.ctypes.data_as(ctypes.c_void_p),
- ao_offsets.ctypes.data_as(ctypes.c_void_p),
- bins_locs_k.ctypes.data_as(ctypes.c_void_p),
- bins_locs_l.ctypes.data_as(ctypes.c_void_p),
- ctypes.c_int(nbins),
- ctypes.c_int(cp_k_id),
- ctypes.c_int(cp_l_id))
-
- if err != 0:
- raise RuntimeError("int2c2e failed\n")
-
- int2c[rows, cols] = int2c[cols, rows]
- if not auxmol.cart:
- coeff = intopt.aux_cart2sph
- int2c = coeff.T @ int2c @ coeff
-
- return int2c
-
-def get_int2c2e(mol, auxmol, direct_scf_tol=1e-13):
- '''
- Generate int2c2e on GPU
- '''
- intopt = VHFOpt(mol, auxmol, 'int2e')
- intopt.build(direct_scf_tol, diag_block_with_triu=True, aosym=True)
- int2c = get_int2c2e_sorted(mol, auxmol, intopt=intopt)
- aux_idx = np.argsort(intopt.aux_ao_idx)
- int2c = int2c[np.ix_(aux_idx, aux_idx)]
- return int2c
-
def sort_mol(mol0, cart=True, log=None):
'''
# Sort basis according to angular momentum and contraction patterns so
diff --git a/gpu4pyscf/lib/gint/nr_fill_ao_ints.cu b/gpu4pyscf/lib/gint/nr_fill_ao_ints.cu
index c77b9814..9d44776b 100644
--- a/gpu4pyscf/lib/gint/nr_fill_ao_ints.cu
+++ b/gpu4pyscf/lib/gint/nr_fill_ao_ints.cu
@@ -142,7 +142,6 @@ int GINTfill_int2e(cudaStream_t stream, BasisProdCache *bpcache, double *eri, in
ContractionProdType *cp_ij = bpcache->cptype + cp_ij_id;
ContractionProdType *cp_kl = bpcache->cptype + cp_kl_id;
GINTEnvVars envs;
-
int ng[4] = {0,0,0,0};
GINTinit_EnvVars(&envs, cp_ij, cp_kl, ng);
envs.omega = omega;
diff --git a/gpu4pyscf/scf/int2c2e.py b/gpu4pyscf/scf/int2c2e.py
new file mode 100644
index 00000000..10b6fc59
--- /dev/null
+++ b/gpu4pyscf/scf/int2c2e.py
@@ -0,0 +1,143 @@
+# Copyright 2024 The GPU4PySCF Authors. All Rights Reserved.
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program. If not, see .
+
+import ctypes
+import cupy
+import numpy as np
+from gpu4pyscf.scf.hf import _make_s_index_offsets
+from gpu4pyscf.lib.cupy_helper import load_library, take_last2d
+from gpu4pyscf.df.int3c2e import VHFOpt, make_fake_mol
+
+libgint = load_library('libgint')
+
+def get_int2c2e_sorted(mol, intopt=None, direct_scf_tol=1e-13, aosym=None, omega=None, stream=None):
+ '''
+ Generated int2c2e consistent with pyscf
+ '''
+ if omega is None: omega = 0.0
+ if stream is None: stream = cupy.cuda.get_current_stream()
+ if intopt is None:
+ # Reuse int3c2e
+ intopt = VHFOpt(mol, mol, 'int2e')
+ intopt.build(direct_scf_tol, diag_block_with_triu=True, aosym=True)
+ nao = mol.nao
+ rows, cols = np.tril_indices(nao)
+
+ nao_cart = intopt.mol.nao
+ norb_cart = nao_cart + 1
+
+ int2c = cupy.zeros([nao_cart, nao_cart], order='F')
+ ao_offsets = np.array([nao_cart+1, nao_cart, nao_cart+1, nao_cart], dtype=np.int32)
+ strides = np.array([1, nao_cart, nao_cart, nao_cart*nao_cart], dtype=np.int32)
+ log_cutoff = np.log(direct_scf_tol)
+ for k_id, log_q_k in enumerate(intopt.aux_log_qs):
+ #bins_locs_k = _make_s_index_offsets(log_q_k, nbins)
+ bins_locs_k = np.array([0, len(log_q_k)], dtype=np.int32)
+ bins_floor_k = np.array([100], dtype=np.double)
+ cp_k_id = k_id + len(intopt.log_qs)
+ for l_id, log_q_l in enumerate(intopt.aux_log_qs):
+ if k_id > l_id: continue
+ #bins_locs_l = _make_s_index_offsets(log_q_l, nbins)
+ bins_locs_l = np.array([0, len(log_q_l)], dtype=np.int32)
+ bins_floor_l = np.array([100], dtype=np.double)
+ cp_l_id = l_id + len(intopt.log_qs)
+ nbins_locs_k = len(bins_locs_k) - 1
+ nbins_locs_l = len(bins_locs_l) - 1
+ err = libgint.GINTfill_int2e(
+ ctypes.cast(stream.ptr, ctypes.c_void_p),
+ intopt.bpcache,
+ ctypes.cast(int2c.data.ptr, ctypes.c_void_p),
+ ctypes.c_int(norb_cart),
+ strides.ctypes.data_as(ctypes.c_void_p),
+ ao_offsets.ctypes.data_as(ctypes.c_void_p),
+ bins_locs_k.ctypes.data_as(ctypes.c_void_p),
+ bins_locs_l.ctypes.data_as(ctypes.c_void_p),
+ bins_floor_k.ctypes.data_as(ctypes.c_void_p),
+ bins_floor_l.ctypes.data_as(ctypes.c_void_p),
+ ctypes.c_int(nbins_locs_k),
+ ctypes.c_int(nbins_locs_l),
+ ctypes.c_int(cp_k_id),
+ ctypes.c_int(cp_l_id),
+ ctypes.c_double(log_cutoff),
+ ctypes.c_double(omega))
+
+ if err != 0:
+ raise RuntimeError("int2c2e failed\n")
+
+ int2c[rows, cols] = int2c[cols, rows]
+ if not mol.cart:
+ coeff = intopt.cart2sph
+ int2c = coeff.T @ int2c @ coeff
+
+ return int2c
+
+def get_int2c2e_ip_sorted(mol, auxmol, intopt=None, direct_scf_tol=1e-13, intor=None, aosym=None, stream=None):
+ '''
+ TODO: WIP
+ '''
+ if stream is None: stream = cupy.cuda.get_current_stream()
+ if intopt is None:
+ intopt = VHFOpt(mol, auxmol, 'int2e')
+ intopt.build(direct_scf_tol, diag_block_with_triu=True, aosym=False)
+
+ nbins = 1
+
+ nao_cart = intopt.mol.nao
+ naux_cart = intopt.auxmol.nao
+ norb_cart = nao_cart + naux_cart + 1
+ rows, cols = np.tril_indices(naux_cart)
+
+ int2c = cupy.zeros([naux_cart, naux_cart], order='F')
+ ao_offsets = np.array([nao_cart+1, nao_cart, nao_cart+1, nao_cart], dtype=np.int32)
+ strides = np.array([1, naux_cart, naux_cart, naux_cart*naux_cart], dtype=np.int32)
+ for k_id, log_q_k in enumerate(intopt.aux_log_qs):
+ bins_locs_k = _make_s_index_offsets(log_q_k, nbins)
+ cp_k_id = k_id + len(intopt.log_qs)
+ for l_id, log_q_l in enumerate(intopt.aux_log_qs):
+ if k_id > l_id: continue
+ bins_locs_l = _make_s_index_offsets(log_q_l, nbins)
+ cp_l_id = l_id + len(intopt.log_qs)
+ err = libgint.GINTfill_int2e(
+ ctypes.cast(stream.ptr, ctypes.c_void_p),
+ intopt.bpcache,
+ ctypes.cast(int2c.data.ptr, ctypes.c_void_p),
+ ctypes.c_int(norb_cart),
+ strides.ctypes.data_as(ctypes.c_void_p),
+ ao_offsets.ctypes.data_as(ctypes.c_void_p),
+ bins_locs_k.ctypes.data_as(ctypes.c_void_p),
+ bins_locs_l.ctypes.data_as(ctypes.c_void_p),
+ ctypes.c_int(nbins),
+ ctypes.c_int(cp_k_id),
+ ctypes.c_int(cp_l_id))
+
+ if err != 0:
+ raise RuntimeError("int2c2e failed\n")
+
+ int2c[rows, cols] = int2c[cols, rows]
+ if not auxmol.cart:
+ coeff = intopt.aux_cart2sph
+ int2c = coeff.T @ int2c @ coeff
+
+ return int2c
+
+def get_int2c2e(mol, direct_scf_tol=1e-13):
+ '''
+ Generate int2c2e on GPU
+ '''
+ intopt = VHFOpt(mol, mol, 'int2e')
+ intopt.build(direct_scf_tol, diag_block_with_triu=True, aosym=True)
+ int2c = get_int2c2e_sorted(mol, intopt=intopt)
+ int2c = take_last2d(int2c, intopt.rev_ao_idx)
+ return int2c
diff --git a/gpu4pyscf/scf/int4c2e.py b/gpu4pyscf/scf/int4c2e.py
index 9636a0ad..8e4049a1 100644
--- a/gpu4pyscf/scf/int4c2e.py
+++ b/gpu4pyscf/scf/int4c2e.py
@@ -438,4 +438,5 @@ def loop_int4c2e_general(intopt, ip_type='', direct_scf_tol=1e-13, omega=None, s
raise RuntimeError(f'GINT_fill_int4c2e general failed, err={err}')
if cp_ij_id == cp_kl_id:
int4c *= 0.5
- yield i0,i1,j0,j1,k0,k1,l0,l1,int4c
\ No newline at end of file
+ yield i0,i1,j0,j1,k0,k1,l0,l1,int4c
+
\ No newline at end of file
diff --git a/gpu4pyscf/scf/tests/test_int2c2e.py b/gpu4pyscf/scf/tests/test_int2c2e.py
new file mode 100644
index 00000000..2a9ae4c5
--- /dev/null
+++ b/gpu4pyscf/scf/tests/test_int2c2e.py
@@ -0,0 +1,76 @@
+# Copyright 2023 The GPU4PySCF Authors. All Rights Reserved.
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program. If not, see .
+
+import pyscf
+from pyscf import lib, df
+from pyscf.scf import _vhf
+from pyscf.gto.moleintor import getints, make_cintopt
+from pyscf.df.grad.rhf import _int3c_wrapper
+import numpy as np
+import cupy
+import unittest
+
+from gpu4pyscf.scf import int2c2e
+from gpu4pyscf.lib.cupy_helper import load_library
+libgint = load_library('libgint')
+
+'''
+compare int2c2e by pyscf and gpu4pyscf
+'''
+
+def setUpModule():
+ global mol_cart, mol_sph
+ mol_cart = pyscf.M(atom='''
+H -0.7570000000 -0.0000000000 -0.4696000000
+H 0.7570000000 0.0000000000 -0.4696000000
+ ''',
+ basis= 'ccpvdz',
+ verbose=1,
+ output = '/dev/null')
+ mol_cart.build()
+
+ mol_sph = pyscf.M(atom='''
+H -0.7570000000 -0.0000000000 -0.4696000000
+H 0.7570000000 0.0000000000 -0.4696000000
+ ''',
+ basis= 'ccpvdz',
+ verbose=1,
+ output = '/dev/null')
+ mol_sph.build()
+
+omega = 0.2
+
+def tearDownModule():
+ global mol_sph, mol_cart
+ mol_sph.stdout.close()
+ mol_cart.stdout.close()
+ del mol_sph, mol_cart
+
+class KnownValues(unittest.TestCase):
+ def test_int2c2e_sph(self):
+ nao = mol_sph.nao
+ eri = mol_sph.intor('int2c2e_sph').reshape((nao,)*2)
+ int2c = int2c2e.get_int2c2e(mol_sph)
+ assert np.linalg.norm(eri - int2c.get()) < 1e-10
+
+ def test_int2c2e_cart(self):
+ nao = mol_cart.nao
+ eri = mol_cart.intor('int2c2e_cart').reshape((nao,)*2)
+ int2c = int2c2e.get_int2c2e(mol_cart)
+ assert np.linalg.norm(eri - int2c.get()) < 1e-10
+
+if __name__ == "__main__":
+ print("Full Tests for int2c")
+ unittest.main()