Assertion failure crash while computing 3-center ERI derivatives
smparker opened this issue · 5 comments
Hi there,
I've run into a crash caused by what looks like an internal assertion failing. I'm still not sure if it's a bug or if I'm just misunderstanding the interface, so I'm hoping you can help figure it out. I'm trying to compute the first derivative of the 3-center ERIs using the c++11 interface. The basic error comes from line 1745 of my engine.impl.h (I realize the line number might change with configure options) and the failing assertion reads buildfnptrs_[buildfnidx] && "null build function ptr".
I have been able to make a minimum example that shows the same behavior that I'll paste below. But first, details.
- I'm using v2.7.1 with a fairly standard configuration (I can share if that'd help)
- It only seems to effect 1st derivatives of 3c integrals.
- regular 3c integrals (no derivative) seem to work fine with
Braket::xs_xx - 2c integrals and their derivatives seem to work fine with
Braket::xs_xs
- regular 3c integrals (no derivative) seem to work fine with
- It doesn't seem to matter whether I use the
engine.compute()orengine.compute2()interfaces - It always crashes at the first p function the code encounters. Increasing
max_ldoesn't seem to change this at all
The test code that triggers this is below.
#include <stdlib.h>
#include <libint2.hpp>
using libint2::BasisSet;
using libint2::Operator;
using libint2::Shell;
int main(int argc, char* argv[])
{
printf("libint2 version %s\n", LIBINT_VERSION);
std::vector<libint2::Shell> shells;
shells.push_back({
{0.100000},
{// S shell, spherical=false, contraction coefficients
{0, false, {1.0}}},
{{0, 0, 0}} // origin coordinates
});
shells.push_back({
{0.100000},
{// P shell, spherical=false, contraction coefficients
{1, false, {1.0}}},
{{0, 0, 0}} // origin coordinates
});
#if LIBINT_MAJOR_VERSION == 2 && LIBINT_MINOR_VERSION >= 8
auto obs = libint2::BasisSet(std::move(shells));
auto shell2bf = obs.shell2bf();
#else
libint2::BasisSet obs = libint2::BasisSet();
for (libint2::Shell shell : shells) {
obs.push_back(std::move(shell));
}
auto shell2bf = BasisSet::compute_shell2bf(obs);
#endif
libint2::initialize();
printf("2-center ERI derivatives\n\n\n");
{ // compute 2-center ERI
libint2::Engine engine = libint2::Engine(libint2::Operator::coulomb, libint2::max_nprim(obs),
libint2::max_l(obs), 1);
engine.set(libint2::BraKet::xs_xs);
const auto& buf = engine.results();
for (auto s1 = 0; s1 != obs.size(); ++s1) {
const auto bf1_first = shell2bf[s1]; // first basis function in this shell
auto n1 = obs[s1].size(); // number of basis functions in shell 1
for (auto s2 = 0; s2 != obs.size(); ++s2) {
const auto bf2_first = shell2bf[s2]; // first basis function in this shell
auto n2 = obs[s2].size(); // number of basis functions in shell 2
engine.compute(obs[s1], obs[s2]);
const auto* buf_1234 = buf[0];
for (auto f1 = 0ul, f1234 = 0ul; f1 != n1; ++f1) {
const auto bf1 = f1 + bf1_first;
for (auto f2 = 0ul; f2 != n2; ++f2, ++f1234) {
const auto bf2 = f2 + bf2_first;
const auto integral = buf_1234[f1234];
printf("d/dx (%lu, %lu), %lf \n", bf1, bf2, integral);
}
}
printf("\n");
}
}
}
printf("3-center ERI derivatives\n\n\n");
{ // compute 3-center ERI
libint2::Engine engine = libint2::Engine(libint2::Operator::coulomb, libint2::max_nprim(obs),
libint2::max_l(obs), 1);
engine.set(libint2::BraKet::xs_xx);
const auto& buf = engine.results();
for (auto s1 = 0; s1 != obs.size(); ++s1) {
const auto bf1_first = shell2bf[s1]; // first basis function in this shell
auto n1 = obs[s1].size(); // number of basis functions in shell 1
for (auto s2 = 0; s2 != obs.size(); ++s2) {
const auto bf2_first = shell2bf[s2]; // first basis function in this shell
auto n2 = obs[s2].size(); // number of basis functions in shell 2
for (auto s3 = 0; s3 != obs.size(); ++s3) {
const auto bf3_first = shell2bf[s3]; // first basis function in this shell
auto n3 = obs[s3].size(); // number of basis functions in shell 3
engine.compute(obs[s1], obs[s2], obs[s3]);
const auto* buf_1234 = buf[0];
for (auto f1 = 0ul, f1234 = 0ul; f1 != n1; ++f1) {
const auto bf1 = f1 + bf1_first;
for (auto f2 = 0ul; f2 != n2; ++f2) {
const auto bf2 = f2 + bf2_first;
for (auto f3 = 0ul; f3 != n3; ++f3, ++f1234) {
const auto bf3 = f3 + bf3_first;
const auto integral = buf_1234[f1234];
printf("d/dx (%lu, %lu, %lu), %lf \n", bf1, bf2, bf3, integral);
}
}
}
printf("\n");
}
}
}
}
libint2::finalize();
return 0;
}
The full output I see from this is:
libint2 version 2.8.1
2-center ERI derivatives
d/dx (0, 0), 0.000000
d/dx (0, 1), 13.246118
d/dx (0, 2), 0.000000
d/dx (0, 3), 0.000000
d/dx (1, 0), -13.246118
d/dx (2, 0), 0.000000
d/dx (3, 0), 0.000000
d/dx (1, 1), 0.000000
d/dx (1, 2), 0.000000
d/dx (1, 3), 0.000000
d/dx (2, 1), 0.000000
d/dx (2, 2), 0.000000
d/dx (2, 3), 0.000000
d/dx (3, 1), 0.000000
d/dx (3, 2), 0.000000
d/dx (3, 3), 0.000000
3-center ERI derivatives
d/dx (0, 0, 0), 0.000000
mintest: /usr/local/include/libint2/engine.impl.h:1875: const target_ptr_vec& libint2::Engine::compute2(const libint2::Shell&, const libint2::Shell&, const libint2::Shell&, const libint2::Shell&, const libint2::ShellPair*, const libint2::ShellPair*) [with libint2::Operator oper = libint2::Operator::coulomb; libint2::BraKet braket = libint2::BraKet::xs_xx; long unsigned int deriv_order = 1; libint2::Engine::target_ptr_vec = std::vector<const double*, libint2::detail::ext_stack_allocator<const double*, 25> >]: Assertion `buildfnptrs_[buildfnidx] && "null build function ptr"' failed.
Aborted
So could this be a bug? Or maybe it is time for me to learn the C interface...
update: I updated the code snippet and results to correspond to 2.8.1
Shane, this is likely a logic error in the Engine ... I will try to reproduce. Unlike 1e 2c and 2e 4c derivatives which have been stressed hard 2e 2c and 2e 3c have not been tested.
So could this be a bug? Or maybe it is time for me to learn the C interface...
NOOOOOOOOOOOOOOOOOOOOOOOO ...
I learned a couple more details about this that I wanted to share just in case they'd be helpful:
- I get basically the same error with 2.8.1 (though I did have to update the example above, which I'll edit shortly)
- It seems that (p|ss) integral derivatives run through, but it crashes at (s|sp) and (s|ps)
Okay, I pretty much know what's going on now. It seems to be caused by using different max angular momentum values for the eri2 and eri3 options. Originally, I used --with-max-am=4 --with-eri-max-am=4 --with-eri2-max-am=6 --with-eri3-max-am=6 when generating the library. When I replace all that with just --with-max-am=6 then the above code seems to work just fine. I can live with that as a solution.