In our implementation, we followed the approach proposed by Ufimtsev and Martinez, where the J and K matrices are formed separately to maximize data locality and to avoid communication between thread blocks. In this approach, the bra and ket pairs are presorted based on the values of Q µν and Q λσ |D λσ |, respectively, as to ensure efficient screening of ERIs in the calculation of the J matrix