Skip to content

refactor: remove Array_Pool in favor of std::vector<double> #7234

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
dzzz2001 wants to merge 2 commits into
base: develop
Choose a base branch
from
Open

refactor: remove Array_Pool in favor of std::vector<double> #7234

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
ddcdf0a
refactor: remove Array_Pool in favor of std::vector<double>
dzzz2001 Apr 10, 2026
4dfdf96
Merge branch 'develop' into remove-array-pool
dzzz2001 Apr 10, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
89 changes: 0 additions & 89 deletions source/source_base/array_pool.h
View file
Open in desktop

This file was deleted.

11 changes: 5 additions & 6 deletions source/source_base/math_ylmreal.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,6 @@
#include "source_base/kernels/math_ylm_op.h"
#include "source_base/libm/libm.h"
#include "source_base/module_device/memory_op.h"
#include "source_base/array_pool.h"
#include "realarray.h"
#include "timer.h"
#include "tool_quit.h"
Expand Down Expand Up @@ -632,7 +631,7 @@ void YlmReal::grad_Ylm_Real
ModuleBase::Ylm::set_coefficients();
const int lmax = int(sqrt( double(lmax2) ) + 0.1) - 1;
std::vector<double> tmpylm((lmax2+1) * (lmax2+1));
Array_Pool<double> tmpgylm((lmax2+1) * (lmax2+1), 3);
std::vector<double> tmpgylm((lmax2+1) * (lmax2+1) * 3);

for (int ig = 0;ig < ng;ig++)
{
Expand All @@ -651,17 +650,17 @@ void YlmReal::grad_Ylm_Real
}
else
{
Ylm::grad_rl_sph_harm(lmax2, gg.x, gg.y, gg.z, tmpylm.data(),tmpgylm.get_ptr_2D());
Ylm::grad_rl_sph_harm(lmax2, gg.x, gg.y, gg.z, tmpylm.data(), tmpgylm.data());
int lm = 0;
for(int il = 0 ; il <= lmax ; ++il)
{
for(int im = 0; im < 2*il+1; ++im, ++lm)
{
double rlylm = tmpylm[lm];
ylm(lm,ig) = rlylm / pow(gmod,il);
dylmx(lm,ig) = ( tmpgylm[lm][0] - il*rlylm * gg.x / pow(gmod,2) )/pow(gmod,il);
dylmy(lm,ig) = ( tmpgylm[lm][1] - il*rlylm * gg.y / pow(gmod,2) )/pow(gmod,il);
dylmz(lm,ig) = ( tmpgylm[lm][2] - il*rlylm * gg.z / pow(gmod,2) )/pow(gmod,il);
dylmx(lm,ig) = ( tmpgylm[lm*3] - il*rlylm * gg.x / pow(gmod,2) )/pow(gmod,il);
dylmy(lm,ig) = ( tmpgylm[lm*3 + 1] - il*rlylm * gg.y / pow(gmod,2) )/pow(gmod,il);
dylmz(lm,ig) = ( tmpgylm[lm*3 + 2] - il*rlylm * gg.z / pow(gmod,2) )/pow(gmod,il);
}
}

Expand Down
15 changes: 7 additions & 8 deletions source/source_base/test/math_ylmreal_test.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,6 @@
#include"gtest/gtest.h"
#include<math.h>
#include "source_psi/psi.h"
#include "source_base/array_pool.h"

#define doublethreshold 1e-12

Expand Down Expand Up @@ -51,7 +50,7 @@ class YlmRealTest : public testing::Test
double *rly; //Ylm
double (*rlgy)[3]; //the gradient of Ylm
std::vector<double> rlyvector; //Ylm
ModuleBase::Array_Pool<double> rlgyvector; //the gradient of Ylm
std::vector<double> rlgyvector; //the gradient of Ylm (flat, size nylm*3)

//Ylm function
inline double norm(const double &x, const double &y, const double &z) {return sqrt(x*x + y*y + z*z);}
Expand Down Expand Up @@ -195,7 +194,7 @@ class YlmRealTest : public testing::Test
rly = new double[nylm];
rlyvector.resize(nylm);
rlgy = new double[nylm][3];
rlgyvector = ModuleBase::Array_Pool<double>(nylm,3);
rlgyvector.assign(nylm * 3, 0.0);
ref = new double[64*4]{
y00(g[0].x, g[0].y, g[0].z), y00(g[1].x, g[1].y, g[1].z), y00(g[2].x, g[2].y, g[2].z), y00(g[3].x, g[3].y, g[3].z),
y10(g[0].x, g[0].y, g[0].z), y10(g[1].x, g[1].y, g[1].z), y10(g[2].x, g[2].y, g[2].z), y10(g[3].x, g[3].y, g[3].z),
Expand Down Expand Up @@ -412,11 +411,11 @@ TEST_F(YlmRealTest,YlmGradRlSphHarm)
for(int j=0;j<ng;++j)
{
double r = sqrt(g[j].x * g[j].x + g[j].y * g[j].y + g[j].z * g[j].z);
ModuleBase::Ylm::grad_rl_sph_harm(lmax,g[j].x/r,g[j].y/r,g[j].z/r,rlyvector.data(),rlgyvector.get_ptr_2D());
ModuleBase::Ylm::grad_rl_sph_harm(lmax,g[j].x/r,g[j].y/r,g[j].z/r,rlyvector.data(),rlgyvector.data());
for(int i=0;i<nylm;++i)
{
EXPECT_NEAR(rlyvector[i],ref[i*ng+j],doublethreshold) << "Ylm[" << i << "], example " << j << " not pass";
for(int k=0;k<3;++k) {EXPECT_NEAR(rlgyvector[i][k],rlgyref[j][i][k],1e-5);}
for(int k=0;k<3;++k) {EXPECT_NEAR(rlgyvector[i*3+k],rlgyref[j][i][k],1e-5);}

}
}
Expand Down Expand Up @@ -467,15 +466,15 @@ TEST_F(YlmRealTest, equality_gradient_test)
double rlya[100];
double rlyb[400];

ModuleBase::Array_Pool<double> grlya(100, 3);
std::vector<double> grlya(100 * 3);
double grlyb[400][3];

ModuleBase::Ylm::grad_rl_sph_harm (9, R.x, R.y, R.z, rlya, grlya.get_ptr_2D());
ModuleBase::Ylm::grad_rl_sph_harm (9, R.x, R.y, R.z, rlya, grlya.data());
ModuleBase::Ylm::rlylm (10, R.x, R.y, R.z, rlyb, grlyb);

for (int i = 0; i < 100; i++)
{
double diffx = fabs(grlya[i][2]-grlyb[i][2]);
double diffx = fabs(grlya[i*3 + 2]-grlyb[i][2]);
EXPECT_LT(diffx,1e-8);
}

Expand Down
138 changes: 38 additions & 100 deletions source/source_base/test/ylm_test.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -77,33 +77,20 @@ TEST_F(ylmTest, HessianFiniteDifferenceL5)
ModuleBase::Ylm::hes_rl_sph_harm(l, x, y, z, hrly);

// Allocate gradient arrays for central difference
std::vector<double> rly_xp((l+1)*(l+1));
std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
double** grly_xp_ptr = new double*[(l+1)*(l+1)];
for (int i = 0; i < (l+1)*(l+1); i++) {
grly_xp_ptr[i] = grly_xp[i].data();
}

std::vector<double> rly_xm((l+1)*(l+1));
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
double** grly_xm_ptr = new double*[(l+1)*(l+1)];
for (int i = 0; i < (l+1)*(l+1); i++) {
grly_xm_ptr[i] = grly_xm[i].data();
}
const int nylm = (l+1)*(l+1);
std::vector<double> rly_xp(nylm), rly_xm(nylm);
std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);

// Compute gradient at (x+h, y, z) and (x-h, y, z)
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm.data());

// Test H_xx for m=0 (index 25) using central difference
int idx = 25;
double H_xx_fd = (grly_xp[idx][0] - grly_xm[idx][0]) / (2.0 * h);
double H_xx_fd = (grly_xp[idx*3] - grly_xm[idx*3]) / (2.0 * h);
double H_xx_analytic = hrly[idx][0];

EXPECT_NEAR(H_xx_fd, H_xx_analytic, tol);

delete[] grly_xp_ptr;
delete[] grly_xm_ptr;
}

// Test Hessian finite difference for l=6 using central difference
Expand All @@ -118,33 +105,20 @@ TEST_F(ylmTest, HessianFiniteDifferenceL6)
ModuleBase::Ylm::hes_rl_sph_harm(l, x, y, z, hrly);

// Allocate gradient arrays for central difference
std::vector<double> rly_xp((l+1)*(l+1));
std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
double** grly_xp_ptr = new double*[(l+1)*(l+1)];
for (int i = 0; i < (l+1)*(l+1); i++) {
grly_xp_ptr[i] = grly_xp[i].data();
}

std::vector<double> rly_xm((l+1)*(l+1));
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
double** grly_xm_ptr = new double*[(l+1)*(l+1)];
for (int i = 0; i < (l+1)*(l+1); i++) {
grly_xm_ptr[i] = grly_xm[i].data();
}
const int nylm = (l+1)*(l+1);
std::vector<double> rly_xp(nylm), rly_xm(nylm);
std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);

// Compute gradient at (x+h, y, z) and (x-h, y, z)
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm.data());

// Test H_xx for m=0 (index 36) using central difference
int idx = 36;
double H_xx_fd = (grly_xp[idx][0] - grly_xm[idx][0]) / (2.0 * h);
double H_xx_fd = (grly_xp[idx*3] - grly_xm[idx*3]) / (2.0 * h);
double H_xx_analytic = hrly[idx][0];

EXPECT_NEAR(H_xx_fd, H_xx_analytic, tol);

delete[] grly_xp_ptr;
delete[] grly_xm_ptr;
}

// Test that l>6 triggers error
Expand Down Expand Up @@ -174,71 +148,46 @@ TEST_F(ylmTest, HessianAllComponentsL2)
int idx = 4;

// Allocate gradient arrays
std::vector<double> rly_xp((l+1)*(l+1)), rly_xm((l+1)*(l+1));
std::vector<double> rly_yp((l+1)*(l+1)), rly_ym((l+1)*(l+1));
std::vector<double> rly_zp((l+1)*(l+1)), rly_zm((l+1)*(l+1));

std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
std::vector<std::vector<double>> grly_yp((l+1)*(l+1), std::vector<double>(3));
std::vector<std::vector<double>> grly_ym((l+1)*(l+1), std::vector<double>(3));
std::vector<std::vector<double>> grly_zp((l+1)*(l+1), std::vector<double>(3));
std::vector<std::vector<double>> grly_zm((l+1)*(l+1), std::vector<double>(3));

double** grly_xp_ptr = new double*[(l+1)*(l+1)];
double** grly_xm_ptr = new double*[(l+1)*(l+1)];
double** grly_yp_ptr = new double*[(l+1)*(l+1)];
double** grly_ym_ptr = new double*[(l+1)*(l+1)];
double** grly_zp_ptr = new double*[(l+1)*(l+1)];
double** grly_zm_ptr = new double*[(l+1)*(l+1)];

for (int i = 0; i < (l+1)*(l+1); i++) {
grly_xp_ptr[i] = grly_xp[i].data();
grly_xm_ptr[i] = grly_xm[i].data();
grly_yp_ptr[i] = grly_yp[i].data();
grly_ym_ptr[i] = grly_ym[i].data();
grly_zp_ptr[i] = grly_zp[i].data();
grly_zm_ptr[i] = grly_zm[i].data();
}
const int nylm = (l+1)*(l+1);
std::vector<double> rly_xp(nylm), rly_xm(nylm);
std::vector<double> rly_yp(nylm), rly_ym(nylm);
std::vector<double> rly_zp(nylm), rly_zm(nylm);

std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);
std::vector<double> grly_yp(nylm * 3), grly_ym(nylm * 3);
std::vector<double> grly_zp(nylm * 3), grly_zm(nylm * 3);

// Compute gradients at perturbed points
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y+h, z, rly_yp.data(), grly_yp_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y-h, z, rly_ym.data(), grly_ym_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z+h, rly_zp.data(), grly_zp_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z-h, rly_zm.data(), grly_zm_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y+h, z, rly_yp.data(), grly_yp.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y-h, z, rly_ym.data(), grly_ym.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z+h, rly_zp.data(), grly_zp.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z-h, rly_zm.data(), grly_zm.data());

// Test H_xx (index 0)
double H_xx_fd = (grly_xp[idx][0] - grly_xm[idx][0]) / (2.0 * h);
double H_xx_fd = (grly_xp[idx*3] - grly_xm[idx*3]) / (2.0 * h);
EXPECT_NEAR(H_xx_fd, hrly[idx][0], tol);

// Test H_xy (index 1)
double H_xy_fd = (grly_xp[idx][1] - grly_xm[idx][1]) / (2.0 * h);
double H_xy_fd = (grly_xp[idx*3 + 1] - grly_xm[idx*3 + 1]) / (2.0 * h);
EXPECT_NEAR(H_xy_fd, hrly[idx][1], tol);

// Test H_xz (index 2)
double H_xz_fd = (grly_xp[idx][2] - grly_xm[idx][2]) / (2.0 * h);
double H_xz_fd = (grly_xp[idx*3 + 2] - grly_xm[idx*3 + 2]) / (2.0 * h);
EXPECT_NEAR(H_xz_fd, hrly[idx][2], tol);

// Test H_yy (index 3)
double H_yy_fd = (grly_yp[idx][1] - grly_ym[idx][1]) / (2.0 * h);
double H_yy_fd = (grly_yp[idx*3 + 1] - grly_ym[idx*3 + 1]) / (2.0 * h);
EXPECT_NEAR(H_yy_fd, hrly[idx][3], tol);

// Test H_yz (index 4)
double H_yz_fd = (grly_yp[idx][2] - grly_ym[idx][2]) / (2.0 * h);
double H_yz_fd = (grly_yp[idx*3 + 2] - grly_ym[idx*3 + 2]) / (2.0 * h);
EXPECT_NEAR(H_yz_fd, hrly[idx][4], tol);

// Test H_zz (index 5)
double H_zz_fd = (grly_zp[idx][2] - grly_zm[idx][2]) / (2.0 * h);
double H_zz_fd = (grly_zp[idx*3 + 2] - grly_zm[idx*3 + 2]) / (2.0 * h);
EXPECT_NEAR(H_zz_fd, hrly[idx][5], tol);

delete[] grly_xp_ptr;
delete[] grly_xm_ptr;
delete[] grly_yp_ptr;
delete[] grly_ym_ptr;
delete[] grly_zp_ptr;
delete[] grly_zm_ptr;
}

// Test Hessian for m=0 values across different l
Expand All @@ -256,28 +205,17 @@ TEST_F(ylmTest, HessianM0DifferentL)
ModuleBase::Ylm::hes_rl_sph_harm(l, x, y, z, hrly);

// Allocate gradient arrays
std::vector<double> rly_xp((l+1)*(l+1)), rly_xm((l+1)*(l+1));
std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
const int nylm = (l+1)*(l+1);
std::vector<double> rly_xp(nylm), rly_xm(nylm);
std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);

double** grly_xp_ptr = new double*[(l+1)*(l+1)];
double** grly_xm_ptr = new double*[(l+1)*(l+1)];

for (int i = 0; i < (l+1)*(l+1); i++) {
grly_xp_ptr[i] = grly_xp[i].data();
grly_xm_ptr[i] = grly_xm[i].data();
}

ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm_ptr);
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp.data());
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm.data());

// Test H_xx for m=0 (index l*l)
int idx = l * l;
double H_xx_fd = (grly_xp[idx][0] - grly_xm[idx][0]) / (2.0 * h);
double H_xx_fd = (grly_xp[idx*3] - grly_xm[idx*3]) / (2.0 * h);
EXPECT_NEAR(H_xx_fd, hrly[idx][0], tol) << "Failed for l=" << l << " m=0";

delete[] grly_xp_ptr;
delete[] grly_xm_ptr;
}
}

Expand Down
Loading
Loading