Skip to content

Commit fbec33a

Browse files
dzzz2001claude
dzzz2001
and
claude
authored
refactor: remove Array_Pool in favor of std::vector<double> (#7234)
Replace ModuleBase::Array_Pool<double> — a bespoke 2D container — with contiguous std::vector<double> throughout the Ylm gradient path. Ylm::grad_rl_sph_harm now takes a flat double* (size n*3) instead of double**. Inside the function, a one-line reinterpret_cast aliases the buffer as double(*)[3] so the 120+ existing grly[lm][xyz] accesses stay untouched — memory layout and codegen are unchanged. All call sites (gint, two-center integrator, math_ylmreal, center2_orb, tests) now allocate std::vector<double>(n*3) and index as grly[lm*3+k]. In set_ddphi.cpp, the constant displ(6,3) table becomes a constexpr double[6][3] literal. array_pool.h is deleted. Co-authored-by: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
1 parent d0a094b commit fbec33a

10 files changed

Lines changed: 90 additions & 241 deletions

File tree

‎source/source_base/array_pool.h‎

Lines changed: 0 additions & 89 deletions
This file was deleted.

‎source/source_base/math_ylmreal.cpp‎

Lines changed: 5 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -4,7 +4,6 @@
44
#include "source_base/kernels/math_ylm_op.h"
55
#include "source_base/libm/libm.h"
66
#include "source_base/module_device/memory_op.h"
7-
#include "source_base/array_pool.h"
87
#include "realarray.h"
98
#include "timer.h"
109
#include "tool_quit.h"
@@ -632,7 +631,7 @@ void YlmReal::grad_Ylm_Real
632631
ModuleBase::Ylm::set_coefficients();
633632
const int lmax = int(sqrt( double(lmax2) ) + 0.1) - 1;
634633
std::vector<double> tmpylm((lmax2+1) * (lmax2+1));
635-
Array_Pool<double> tmpgylm((lmax2+1) * (lmax2+1), 3);
634+
std::vector<double> tmpgylm((lmax2+1) * (lmax2+1) * 3);
636635

637636
for (int ig = 0;ig < ng;ig++)
638637
{
@@ -651,17 +650,17 @@ void YlmReal::grad_Ylm_Real
651650
}
652651
else
653652
{
654-
Ylm::grad_rl_sph_harm(lmax2, gg.x, gg.y, gg.z, tmpylm.data(),tmpgylm.get_ptr_2D());
653+
Ylm::grad_rl_sph_harm(lmax2, gg.x, gg.y, gg.z, tmpylm.data(), tmpgylm.data());
655654
int lm = 0;
656655
for(int il = 0 ; il <= lmax ; ++il)
657656
{
658657
for(int im = 0; im < 2*il+1; ++im, ++lm)
659658
{
660659
double rlylm = tmpylm[lm];
661660
ylm(lm,ig) = rlylm / pow(gmod,il);
662-
dylmx(lm,ig) = ( tmpgylm[lm][0] - il*rlylm * gg.x / pow(gmod,2) )/pow(gmod,il);
663-
dylmy(lm,ig) = ( tmpgylm[lm][1] - il*rlylm * gg.y / pow(gmod,2) )/pow(gmod,il);
664-
dylmz(lm,ig) = ( tmpgylm[lm][2] - il*rlylm * gg.z / pow(gmod,2) )/pow(gmod,il);
661+
dylmx(lm,ig) = ( tmpgylm[lm*3] - il*rlylm * gg.x / pow(gmod,2) )/pow(gmod,il);
662+
dylmy(lm,ig) = ( tmpgylm[lm*3 + 1] - il*rlylm * gg.y / pow(gmod,2) )/pow(gmod,il);
663+
dylmz(lm,ig) = ( tmpgylm[lm*3 + 2] - il*rlylm * gg.z / pow(gmod,2) )/pow(gmod,il);
665664
}
666665
}
667666

‎source/source_base/test/math_ylmreal_test.cpp‎

Lines changed: 7 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -5,7 +5,6 @@
55
#include"gtest/gtest.h"
66
#include<math.h>
77
#include "source_psi/psi.h"
8-
#include "source_base/array_pool.h"
98

109
#define doublethreshold 1e-12
1110

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

5655
//Ylm function
5756
inline double norm(const double &x, const double &y, const double &z) {return sqrt(x*x + y*y + z*z);}
@@ -195,7 +194,7 @@ class YlmRealTest : public testing::Test
195194
rly = new double[nylm];
196195
rlyvector.resize(nylm);
197196
rlgy = new double[nylm][3];
198-
rlgyvector = ModuleBase::Array_Pool<double>(nylm,3);
197+
rlgyvector.assign(nylm * 3, 0.0);
199198
ref = new double[64*4]{
200199
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),
201200
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),
@@ -412,11 +411,11 @@ TEST_F(YlmRealTest,YlmGradRlSphHarm)
412411
for(int j=0;j<ng;++j)
413412
{
414413
double r = sqrt(g[j].x * g[j].x + g[j].y * g[j].y + g[j].z * g[j].z);
415-
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());
414+
ModuleBase::Ylm::grad_rl_sph_harm(lmax,g[j].x/r,g[j].y/r,g[j].z/r,rlyvector.data(),rlgyvector.data());
416415
for(int i=0;i<nylm;++i)
417416
{
418417
EXPECT_NEAR(rlyvector[i],ref[i*ng+j],doublethreshold) << "Ylm[" << i << "], example " << j << " not pass";
419-
for(int k=0;k<3;++k) {EXPECT_NEAR(rlgyvector[i][k],rlgyref[j][i][k],1e-5);}
418+
for(int k=0;k<3;++k) {EXPECT_NEAR(rlgyvector[i*3+k],rlgyref[j][i][k],1e-5);}
420419

421420
}
422421
}
@@ -467,15 +466,15 @@ TEST_F(YlmRealTest, equality_gradient_test)
467466
double rlya[100];
468467
double rlyb[400];
469468

470-
ModuleBase::Array_Pool<double> grlya(100, 3);
469+
std::vector<double> grlya(100 * 3);
471470
double grlyb[400][3];
472471

473-
ModuleBase::Ylm::grad_rl_sph_harm (9, R.x, R.y, R.z, rlya, grlya.get_ptr_2D());
472+
ModuleBase::Ylm::grad_rl_sph_harm (9, R.x, R.y, R.z, rlya, grlya.data());
474473
ModuleBase::Ylm::rlylm (10, R.x, R.y, R.z, rlyb, grlyb);
475474

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

‎source/source_base/test/ylm_test.cpp‎

Lines changed: 38 additions & 100 deletions
Original file line numberDiff line numberDiff line change
@@ -77,33 +77,20 @@ TEST_F(ylmTest, HessianFiniteDifferenceL5)
7777
ModuleBase::Ylm::hes_rl_sph_harm(l, x, y, z, hrly);
7878

7979
// Allocate gradient arrays for central difference
80-
std::vector<double> rly_xp((l+1)*(l+1));
81-
std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
82-
double** grly_xp_ptr = new double*[(l+1)*(l+1)];
83-
for (int i = 0; i < (l+1)*(l+1); i++) {
84-
grly_xp_ptr[i] = grly_xp[i].data();
85-
}
86-
87-
std::vector<double> rly_xm((l+1)*(l+1));
88-
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
89-
double** grly_xm_ptr = new double*[(l+1)*(l+1)];
90-
for (int i = 0; i < (l+1)*(l+1); i++) {
91-
grly_xm_ptr[i] = grly_xm[i].data();
92-
}
80+
const int nylm = (l+1)*(l+1);
81+
std::vector<double> rly_xp(nylm), rly_xm(nylm);
82+
std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);
9383

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

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

10393
EXPECT_NEAR(H_xx_fd, H_xx_analytic, tol);
104-
105-
delete[] grly_xp_ptr;
106-
delete[] grly_xm_ptr;
10794
}
10895

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

120107
// Allocate gradient arrays for central difference
121-
std::vector<double> rly_xp((l+1)*(l+1));
122-
std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
123-
double** grly_xp_ptr = new double*[(l+1)*(l+1)];
124-
for (int i = 0; i < (l+1)*(l+1); i++) {
125-
grly_xp_ptr[i] = grly_xp[i].data();
126-
}
127-
128-
std::vector<double> rly_xm((l+1)*(l+1));
129-
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
130-
double** grly_xm_ptr = new double*[(l+1)*(l+1)];
131-
for (int i = 0; i < (l+1)*(l+1); i++) {
132-
grly_xm_ptr[i] = grly_xm[i].data();
133-
}
108+
const int nylm = (l+1)*(l+1);
109+
std::vector<double> rly_xp(nylm), rly_xm(nylm);
110+
std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);
134111

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

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

144121
EXPECT_NEAR(H_xx_fd, H_xx_analytic, tol);
145-
146-
delete[] grly_xp_ptr;
147-
delete[] grly_xm_ptr;
148122
}
149123

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

176150
// Allocate gradient arrays
177-
std::vector<double> rly_xp((l+1)*(l+1)), rly_xm((l+1)*(l+1));
178-
std::vector<double> rly_yp((l+1)*(l+1)), rly_ym((l+1)*(l+1));
179-
std::vector<double> rly_zp((l+1)*(l+1)), rly_zm((l+1)*(l+1));
180-
181-
std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
182-
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
183-
std::vector<std::vector<double>> grly_yp((l+1)*(l+1), std::vector<double>(3));
184-
std::vector<std::vector<double>> grly_ym((l+1)*(l+1), std::vector<double>(3));
185-
std::vector<std::vector<double>> grly_zp((l+1)*(l+1), std::vector<double>(3));
186-
std::vector<std::vector<double>> grly_zm((l+1)*(l+1), std::vector<double>(3));
187-
188-
double** grly_xp_ptr = new double*[(l+1)*(l+1)];
189-
double** grly_xm_ptr = new double*[(l+1)*(l+1)];
190-
double** grly_yp_ptr = new double*[(l+1)*(l+1)];
191-
double** grly_ym_ptr = new double*[(l+1)*(l+1)];
192-
double** grly_zp_ptr = new double*[(l+1)*(l+1)];
193-
double** grly_zm_ptr = new double*[(l+1)*(l+1)];
194-
195-
for (int i = 0; i < (l+1)*(l+1); i++) {
196-
grly_xp_ptr[i] = grly_xp[i].data();
197-
grly_xm_ptr[i] = grly_xm[i].data();
198-
grly_yp_ptr[i] = grly_yp[i].data();
199-
grly_ym_ptr[i] = grly_ym[i].data();
200-
grly_zp_ptr[i] = grly_zp[i].data();
201-
grly_zm_ptr[i] = grly_zm[i].data();
202-
}
151+
const int nylm = (l+1)*(l+1);
152+
std::vector<double> rly_xp(nylm), rly_xm(nylm);
153+
std::vector<double> rly_yp(nylm), rly_ym(nylm);
154+
std::vector<double> rly_zp(nylm), rly_zm(nylm);
155+
156+
std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);
157+
std::vector<double> grly_yp(nylm * 3), grly_ym(nylm * 3);
158+
std::vector<double> grly_zp(nylm * 3), grly_zm(nylm * 3);
203159

204160
// Compute gradients at perturbed points
205-
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp_ptr);
206-
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm_ptr);
207-
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y+h, z, rly_yp.data(), grly_yp_ptr);
208-
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y-h, z, rly_ym.data(), grly_ym_ptr);
209-
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z+h, rly_zp.data(), grly_zp_ptr);
210-
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z-h, rly_zm.data(), grly_zm_ptr);
161+
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp.data());
162+
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm.data());
163+
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y+h, z, rly_yp.data(), grly_yp.data());
164+
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y-h, z, rly_ym.data(), grly_ym.data());
165+
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z+h, rly_zp.data(), grly_zp.data());
166+
ModuleBase::Ylm::grad_rl_sph_harm(l, x, y, z-h, rly_zm.data(), grly_zm.data());
211167

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

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

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

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

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

232188
// Test H_zz (index 5)
233-
double H_zz_fd = (grly_zp[idx][2] - grly_zm[idx][2]) / (2.0 * h);
189+
double H_zz_fd = (grly_zp[idx*3 + 2] - grly_zm[idx*3 + 2]) / (2.0 * h);
234190
EXPECT_NEAR(H_zz_fd, hrly[idx][5], tol);
235-
236-
delete[] grly_xp_ptr;
237-
delete[] grly_xm_ptr;
238-
delete[] grly_yp_ptr;
239-
delete[] grly_ym_ptr;
240-
delete[] grly_zp_ptr;
241-
delete[] grly_zm_ptr;
242191
}
243192

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

258207
// Allocate gradient arrays
259-
std::vector<double> rly_xp((l+1)*(l+1)), rly_xm((l+1)*(l+1));
260-
std::vector<std::vector<double>> grly_xp((l+1)*(l+1), std::vector<double>(3));
261-
std::vector<std::vector<double>> grly_xm((l+1)*(l+1), std::vector<double>(3));
208+
const int nylm = (l+1)*(l+1);
209+
std::vector<double> rly_xp(nylm), rly_xm(nylm);
210+
std::vector<double> grly_xp(nylm * 3), grly_xm(nylm * 3);
262211

263-
double** grly_xp_ptr = new double*[(l+1)*(l+1)];
264-
double** grly_xm_ptr = new double*[(l+1)*(l+1)];
265-
266-
for (int i = 0; i < (l+1)*(l+1); i++) {
267-
grly_xp_ptr[i] = grly_xp[i].data();
268-
grly_xm_ptr[i] = grly_xm[i].data();
269-
}
270-
271-
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp_ptr);
272-
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm_ptr);
212+
ModuleBase::Ylm::grad_rl_sph_harm(l, x+h, y, z, rly_xp.data(), grly_xp.data());
213+
ModuleBase::Ylm::grad_rl_sph_harm(l, x-h, y, z, rly_xm.data(), grly_xm.data());
273214

274215
// Test H_xx for m=0 (index l*l)
275216
int idx = l * l;
276-
double H_xx_fd = (grly_xp[idx][0] - grly_xm[idx][0]) / (2.0 * h);
217+
double H_xx_fd = (grly_xp[idx*3] - grly_xm[idx*3]) / (2.0 * h);
277218
EXPECT_NEAR(H_xx_fd, hrly[idx][0], tol) << "Failed for l=" << l << " m=0";
278-
279-
delete[] grly_xp_ptr;
280-
delete[] grly_xm_ptr;
281219
}
282220
}
283221

0 commit comments

Comments
 (0)