fuse_auxcell assumes l ≤ 7 and crashes with high-angular-momentum auxiliary basis

[WSLinkK]

Description

When using green_igen with an auxiliary cell that contains high-angular-momentum functions (e.g. g functions), the code crashes inside fuse_auxcell with an IndexError. The function currently assumes that the maximum angular momentum in the compensating charge basis (chgcell) is at most 7, and allocates an array with a hard-coded second dimension of size 8. Same issue has been noted by Pavel in upstream in PySCF in issues [pyscf/pyscf#22169]

Relevant snippet from green_igen/df.py:

def fuse_auxcell(mydf, auxcell):
    eta = mydf.eta
    # if eta is None:
    #     eta, mesh, ke_cutoff = _guess_eta(auxcell, mydf.kpts, mydf.mesh)
    chgcell = make_modchg_basis(auxcell, eta)
    fused_cell = copy.copy(auxcell)
    fused_cell._atm, fused_cell._bas, fused_cell._env = \
            gto.conc_env(auxcell._atm, auxcell._bas, auxcell._env,
                         chgcell._atm, chgcell._bas, chgcell._env)
    fused_cell.rcut = max(auxcell.rcut, chgcell.rcut)

    aux_loc = auxcell.ao_loc_nr()
    naux = aux_loc[-1]
    modchg_offset = -numpy.ones((chgcell.natm, 8), dtype=int)  # <-- hard-coded 8
    smooth_loc = chgcell.ao_loc_nr()
    for i in range(chgcell.nbas):
        ia = chgcell.bas_atom(i)
        l  = chgcell.bas_angular(i)
        modchg_offset[ia, l] = smooth_loc[i]

Error message

Traceback (most recent call last):
  File "/home/munkhw/lib/green/python/init_data_df.py", line 9, in <module>
    pyscf_init.mean_field_input()
  File "/home/munkhw/lib/green-mbtools/green_mbtools/mint/pyscf_init.py", line 137, in mean_field_input
    self.compute_df_int(nao, X_k)
  File "/home/munkhw/lib/green-mbtools/green_mbtools/mint/pyscf_init.py", line 170, in compute_df_int
    gggdf._make_j3c(mydf, self.cell, auxcell, kptij_lst, "cderi_ewald.h5")
  File "/home/munkhw/.local/lib/python3.10/site-packages/green_igen/df.py", line 219, in _make_j3c
    fused_cell, fuse = fuse_auxcell(mydf, auxcell)
  File "/home/munkhw/.local/lib/python3.10/site-packages/green_igen/df.py", line 1129, in fuse_auxcell
    modchg_offset[ia,l] = smooth_loc[i]
IndexError: index 8 is out of bounds for axis 1 with size 8

Steps to reproduce

1. Set up a periodic DF calculation using gggdf._make_j3c with an auxiliary / main basis that includes g functions (high angular momentum).
2. Run the calculation so that _make_j3c calls fuse_auxcell(mydf, auxcell).
3. The code crashes in the modchg_offset[ia, l] assignment when bas_angular(i) returns a value beyond the hard-coded limit.

Expected behavior

fuse_auxcell should work for arbitrary angular momentum in chgcell and not assume l ≤ 7. The modchg_offset array should be sized according to the actual maximum angular momentum present in the (auxiliary / compensating) cell.

How PySCF fixes it

The same problem was fixed upstream in PySCF in PR [pyscf/pyscf#2361](pyscf/pyscf#2361), in pyscf/pbc/df/gdf_builder.py. The fix replaces the hard-coded dimension 8 by computing the maximum angular momentum and sizing the offset array accordingly, e.g.

aux_loc = auxcell.ao_loc_nr()
naux = aux_loc[-1]
lmax = auxcell._bas[:, gto.ANG_OF].max()
modchg_offset = -np.ones((chgcell.natm, lmax + 1), dtype=int)
smooth_loc = chgcell.ao_loc_nr()

That is:

1. Compute lmax from the basis.
2. Allocate modchg_offset with shape (natm, lmax + 1) so all angular momenta up to the maximum are covered.

Environment

• Python: 3.10.10
• green_igen version: 0.2.7
• PySCF version: 2.7.0
• System / cluster: pauli