Skip to content

[Codescan] Release FLINT/ARB temporaries in every 2F2 evaluation #31

Description

@Stardust0831

Global scan provenance

This issue was found by the Codex global repository scan of every non-unit-test file in Stardust0831/Multiwfn, pinned to commit d2b770a531eb7039153198d04520b5992f15c49d.

Problem

hyp2F2_prec() allocates two ACB vectors and initializes ARB, ACB and ARF objects on every call, but returns without clearing any of them. The early non-real return leaks the same objects.

Multiwfn/ext/2F2.c

Lines 12 to 39 in d2b770a

double hyp2F2_prec(double a1, double a2, double b1, double b2, double z,
slong prec) {
const slong p = 2, q = 2;
const int regularized = 0;
arb_t res_r;
acb_t res_ri, _z;
acb_ptr a, b;
arf_t resf;
a = _acb_vec_init(p);
b = _acb_vec_init(q);
arb_init(res_r);
acb_init(res_ri);
acb_init(_z);
arf_init(resf);
acb_set_d(a + 0, a1);
acb_set_d(a + 1, a2);
acb_set_d(b + 0, b1);
acb_set_d(b + 1, b2);
acb_set_d(_z, z);
acb_hypgeom_pfq(res_ri, a, p, b, q, _z, regularized, prec);
if (acb_is_real(res_ri)) {
acb_get_real(res_r, res_ri);
} else {
return 0.0;
}
arb_get_abs_ubound_arf(resf, res_r, prec);
double hyp2F2 = arf_get_d(resf, ARF_RND_DOWN);
return hyp2F2;

hyp2F2() may call hyp2F2_prec() repeatedly while increasing precision, multiplying the leak for difficult values:

Multiwfn/ext/2F2.c

Lines 42 to 51 in d2b770a

double hyp2F2(double a1, double a2, double b1, double b2, double z) {
slong prec = 1000;
double hyp2F2 = 0.0;
// for our purpose, 0 < hyp2F2 <= 1
while (hyp2F2 == 0.0 || hyp2F2 > 1.0) {
hyp2F2 = hyp2F2_prec(a1, a2, b1, b2, z, prec);
if(hyp2F2 == 0.0 && prec > 100000) break;
prec *= 10;
}
return hyp2F2;

This is not a one-time startup allocation. Fractional orbital evaluation calls GTO_fractional_integral() inside the loop over primitive Gaussian functions, and that path reaches hyp2F2() for general exact evaluation:

Multiwfn/function.f90

Lines 1690 to 1720 in d2b770a

! main cycle over shells
do j = 1, nprims_t
ix = type2ix(b_t(j)%type)
iy = type2iy(b_t(j)%type)
iz = type2iz(b_t(j)%type)
it = ix + iy + iz
ep = b_t(j)%exp
if (b_t(j)%center /= last_center) then
last_center = b_t(j)%center
sftx1 = x - (a(b_t(j)%center)%x + xmove)
sfty1 = y - (a(b_t(j)%center)%y + ymove)
sftz1 = z - (a(b_t(j)%center)%z + zmove)
sftx2 = sftx1 * sftx1
sfty2 = sfty1 * sfty1
sftz2 = sftz1 * sftz1
rr = sftx2 + sfty2 + sftz2
end if
! commented out since cutoff should depend on p-alpha
! Fox example, when p-alpha = 0, expcutoff can be -40 and results are fine
! but when p-alpha = 2, expcutoff should be more than 20000 for avoiding numerical problems
! if (expcutoff>0._8 .or. -ep*rr>-20000._8) then
if (expcutoff>0._8 .or. .true.) then ! all orbitals are important
else
cycle ! just skip it
end if
do ip = 0, ipmax
t(ip) = GTO_fractional_integral(it + ip * 2, -ep*rr)
end do

Multiwfn/ext/2F2.f90

Lines 131 to 137 in d2b770a

if (algorithm_ == GENERAL_EXACT_EVALUATION) then
top = real(n + 1, kind=8)
down = top + real(p, kind=8) - alpha_
gtop = gamma(top)
gdown = gamma(down)
val2F2 = hyp2F2(top, top + one, down, down + one, x)
GTO_fractional_integral = gtop / gdown * val2F2

Impact

Builds configured with MULTIWFN_WITH_FD=ON leak GMP/FLINT-managed memory for every evaluated primitive and potentially for every precision retry. Grid or repeated point evaluation can therefore grow memory continuously and eventually become slow or terminate from memory exhaustion.

The default CI configuration has fractional derivatives disabled, so existing CI does not exercise or detect this path.

Suggested direction

  • Use a single cleanup path before every return.
  • Call _acb_vec_clear(a, p) and _acb_vec_clear(b, q).
  • Call arb_clear(res_r), acb_clear(res_ri), acb_clear(_z) and arf_clear(resf).
  • Preserve the computed double in a plain local, clean all FLINT objects, and only then return.
  • Add an FD-enabled stress/regression test that evaluates many points and checks that resident memory does not grow proportionally to the number of evaluations.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions