|
| 1 | +#include "gtest/gtest.h" |
| 2 | +#include "gmock/gmock.h" |
| 3 | +#include "../module_ewald/H_Ewald_pw.h" |
| 4 | +#include "../module_ewald/dnrm2.h" |
| 5 | +#include "source_base/matrix3.h" |
| 6 | +#include <vector> |
| 7 | + |
| 8 | +/************************************************ |
| 9 | + * unit test of H_Ewald_pw::rgen |
| 10 | + ***********************************************/ |
| 11 | + |
| 12 | +/** |
| 13 | + * - Tested Functions: |
| 14 | + * - H_Ewald_pw::rgen(): |
| 15 | + * - Generates lattice vectors R such that |R - dtau| <= rmax, |
| 16 | + * and returns them sorted by ascending magnitude. |
| 17 | + * - Tested cases: |
| 18 | + * 1. rmax == 0.0: no vectors returned. |
| 19 | + * 2. Simple cubic cell, small rmax: correct count + sorted order. |
| 20 | + * 3. Large rmax exceeding original fixed mxr=200 limit: verifies |
| 21 | + * that the dynamic mxr sizing introduced in the bug fix works |
| 22 | + * correctly and does not overflow the allocated arrays. |
| 23 | + */ |
| 24 | + |
| 25 | +class RgenTest : public ::testing::Test |
| 26 | +{ |
| 27 | + protected: |
| 28 | + // Simple cubic unit cell: latvec = G = identity |
| 29 | + ModuleBase::Matrix3 latvec; |
| 30 | + ModuleBase::Matrix3 G; |
| 31 | + ModuleBase::Vector3<double> dtau; |
| 32 | + |
| 33 | + void SetUp() override |
| 34 | + { |
| 35 | + // Unit cubic cell: direct and reciprocal lattice vectors are identity |
| 36 | + latvec = ModuleBase::Matrix3(1, 0, 0, 0, 1, 0, 0, 0, 1); |
| 37 | + G = ModuleBase::Matrix3(1, 0, 0, 0, 1, 0, 0, 0, 1); |
| 38 | + dtau = ModuleBase::Vector3<double>(0.0, 0.0, 0.0); |
| 39 | + } |
| 40 | +}; |
| 41 | + |
| 42 | +TEST_F(RgenTest, ZeroRmax) |
| 43 | +{ |
| 44 | + // When rmax==0 the function should return immediately with nrm=0 |
| 45 | + const int mxr_test = 10; |
| 46 | + H_Ewald_pw::mxr = mxr_test; |
| 47 | + std::vector<ModuleBase::Vector3<double>> r(mxr_test); |
| 48 | + std::vector<double> r2(mxr_test); |
| 49 | + std::vector<int> irr(mxr_test); |
| 50 | + int nrm = 0; |
| 51 | + |
| 52 | + H_Ewald_pw::rgen(dtau, 0.0, irr.data(), latvec, G, r.data(), r2.data(), nrm); |
| 53 | + |
| 54 | + EXPECT_EQ(nrm, 0); |
| 55 | +} |
| 56 | + |
| 57 | +TEST_F(RgenTest, SimpleCubicNearestNeighbors) |
| 58 | +{ |
| 59 | + // rmax = 1.5 captures nearest (d=1) and next-nearest (d=sqrt(2)~1.414) |
| 60 | + // neighbors: 6 + 12 = 18 vectors total. |
| 61 | + const double rmax = 1.5; |
| 62 | + const int mxr_test = 50; |
| 63 | + H_Ewald_pw::mxr = mxr_test; |
| 64 | + std::vector<ModuleBase::Vector3<double>> r(mxr_test); |
| 65 | + std::vector<double> r2(mxr_test); |
| 66 | + std::vector<int> irr(mxr_test); |
| 67 | + int nrm = 0; |
| 68 | + |
| 69 | + H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), nrm); |
| 70 | + |
| 71 | + EXPECT_EQ(nrm, 18); |
| 72 | + |
| 73 | + // Vectors must be sorted in ascending order of |r|^2 |
| 74 | + for (int i = 1; i < nrm; ++i) |
| 75 | + { |
| 76 | + EXPECT_LE(r2[i - 1], r2[i]); |
| 77 | + } |
| 78 | + |
| 79 | + // All returned vectors must lie strictly inside the sphere |
| 80 | + for (int i = 0; i < nrm; ++i) |
| 81 | + { |
| 82 | + EXPECT_LE(r2[i], rmax * rmax + 1.0e-10); |
| 83 | + EXPECT_GT(r2[i], 1.0e-10); |
| 84 | + } |
| 85 | +} |
| 86 | + |
| 87 | +TEST_F(RgenTest, SimpleCubicNonZeroDtau) |
| 88 | +{ |
| 89 | + // rgen computes t = R - dtau for each lattice vector R=(i,j,k)*latvec, |
| 90 | + // and excludes vectors with |t|^2 < 1e-10 (i.e. R == dtau). |
| 91 | + // With dtau=(0.5,0,0) and rmax=0.6, two vectors qualify: |
| 92 | + // R=(0,0,0): t = (0,0,0)-(0.5,0,0) = (-0.5,0,0), |t|^2=0.25 <= 0.36 |
| 93 | + // R=(1,0,0): t = (1,0,0)-(0.5,0,0) = ( 0.5,0,0), |t|^2=0.25 <= 0.36 |
| 94 | + // No lattice point coincides with dtau, so neither is excluded. |
| 95 | + const double rmax = 0.6; |
| 96 | + const int mxr_test = 10; |
| 97 | + H_Ewald_pw::mxr = mxr_test; |
| 98 | + dtau = ModuleBase::Vector3<double>(0.5, 0.0, 0.0); |
| 99 | + std::vector<ModuleBase::Vector3<double>> r(mxr_test); |
| 100 | + std::vector<double> r2(mxr_test); |
| 101 | + std::vector<int> irr(mxr_test); |
| 102 | + int nrm = 0; |
| 103 | + |
| 104 | + H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), nrm); |
| 105 | + |
| 106 | + EXPECT_EQ(nrm, 2); |
| 107 | + for (int i = 0; i < nrm; ++i) |
| 108 | + { |
| 109 | + EXPECT_NEAR(r2[i], 0.25, 1.0e-10); |
| 110 | + } |
| 111 | +} |
| 112 | + |
| 113 | +TEST_F(RgenTest, LargeRmaxExceedsOriginalLimit) |
| 114 | +{ |
| 115 | + // rmax=4.0 on a unit cubic cell yields ~499 r-vectors, well above the |
| 116 | + // old fixed limit of mxr=200 that caused the buffer overflow. |
| 117 | + // This test verifies that with a properly sized mxr the function |
| 118 | + // completes without error. |
| 119 | + const double rmax = 4.0; |
| 120 | + |
| 121 | + // Replicate the dynamic mxr computation from compute_ewald() |
| 122 | + double bg1[3]; |
| 123 | + bg1[0] = G.e11; bg1[1] = G.e12; bg1[2] = G.e13; |
| 124 | + int nm1 = (int)(dnrm2(3, bg1, 1) * rmax + 2); |
| 125 | + bg1[0] = G.e21; bg1[1] = G.e22; bg1[2] = G.e23; |
| 126 | + int nm2 = (int)(dnrm2(3, bg1, 1) * rmax + 2); |
| 127 | + bg1[0] = G.e31; bg1[1] = G.e32; bg1[2] = G.e33; |
| 128 | + int nm3 = (int)(dnrm2(3, bg1, 1) * rmax + 2); |
| 129 | + const int mxr_test = (2 * nm1 + 1) * (2 * nm2 + 1) * (2 * nm3 + 1); |
| 130 | + |
| 131 | + H_Ewald_pw::mxr = mxr_test; |
| 132 | + std::vector<ModuleBase::Vector3<double>> r(mxr_test); |
| 133 | + std::vector<double> r2(mxr_test); |
| 134 | + std::vector<int> irr(mxr_test); |
| 135 | + int nrm = 0; |
| 136 | + |
| 137 | + H_Ewald_pw::rgen(dtau, rmax, irr.data(), latvec, G, r.data(), r2.data(), nrm); |
| 138 | + |
| 139 | + // Must exceed the old hard-coded limit that caused the crash |
| 140 | + EXPECT_GT(nrm, 200); |
| 141 | + |
| 142 | + // All returned vectors lie within the sphere |
| 143 | + for (int i = 0; i < nrm; ++i) |
| 144 | + { |
| 145 | + EXPECT_LE(r2[i], rmax * rmax + 1.0e-10); |
| 146 | + EXPECT_GT(r2[i], 1.0e-10); |
| 147 | + } |
| 148 | + |
| 149 | + // Vectors are sorted in ascending order of |r|^2 |
| 150 | + for (int i = 1; i < nrm; ++i) |
| 151 | + { |
| 152 | + EXPECT_LE(r2[i - 1], r2[i]); |
| 153 | + } |
| 154 | +} |
0 commit comments