|
| 1 | +#include"../diago_david.h" |
| 2 | +#include"../diago_cg.h" |
| 3 | +#include"./diago_mock.h" |
| 4 | +#include"gtest/gtest.h" |
| 5 | +#include"mpi.h" |
| 6 | + |
| 7 | +#include "../../module_base/inverse_matrix.h" |
| 8 | +#include "../../module_base/lapack_connector.h" |
| 9 | + |
| 10 | +#define CONVTHRESHOLD 1e-3 |
| 11 | +#define DETAILINFO false |
| 12 | + |
| 13 | + |
| 14 | +/************************************************ |
| 15 | +* unit test of class Diago_David |
| 16 | +***********************************************/ |
| 17 | + |
| 18 | +/** |
| 19 | + * Class Diago_David is used to solve the eigenvalues |
| 20 | + * This unittest test the function Diago_David::diag() |
| 21 | + * with different examples. |
| 22 | + * - the hamilt matrix (npw=100,500,1000) produced by random with sparsity of 90% |
| 23 | + * - the hamilt matrix (npw=100,500,1000) produced by random with sparsity of 50% |
| 24 | + * - the hamilt matrix (npw=100,500,1000) produced by random with sparsity of 0% |
| 25 | + * - the hamilt matrix read from "data-H" |
| 26 | + * |
| 27 | + * The test is passed when the eignvalues are closed to these calculated by LAPACK. |
| 28 | + * |
| 29 | + */ |
| 30 | + |
| 31 | +//mock the ddot_real function |
| 32 | +double Diago_CG::ddot_real(int const&, std::complex<double> const*, std::complex<double> const*, bool) {}; |
| 33 | + |
| 34 | +//use lapack to calcualte eigenvalue of matrix hm |
| 35 | +//NOTE: after finish this function, hm stores the eigen vectors. |
| 36 | +void lapackEigen(int &npw, ModuleBase::ComplexMatrix &hm, ModuleBase::ComplexMatrix &ev,double * e, bool outtime=false) |
| 37 | +{ |
| 38 | + int lwork = 2 * npw; |
| 39 | + std::complex<double> *work2= new std::complex<double>[lwork]; |
| 40 | + double* rwork = new double[3*npw-2]; |
| 41 | + int info = 0; |
| 42 | + |
| 43 | + ModuleBase::ComplexMatrix tmp = hm; |
| 44 | + |
| 45 | + clock_t start,end; |
| 46 | + start = clock(); |
| 47 | + LapackConnector::zheev('V', 'U', npw, tmp, npw, e, work2, lwork, rwork, &info); |
| 48 | + end = clock(); |
| 49 | + if (outtime) std::cout<<"Lapack Run time: "<<(double)(end - start) / CLOCKS_PER_SEC<<" S"<<std::endl; |
| 50 | + |
| 51 | + ev = tmp; |
| 52 | + |
| 53 | + delete [] rwork; |
| 54 | + delete [] work2; |
| 55 | +} |
| 56 | + |
| 57 | +class DiagoDavPrepare |
| 58 | +{ |
| 59 | +public: |
| 60 | + DiagoDavPrepare(int nband, int npw, int sparsity, int order,double eps,int maxiter): |
| 61 | + nband(nband),npw(npw),sparsity(sparsity),order(order),eps(eps),maxiter(maxiter) {} |
| 62 | + |
| 63 | + int nband, npw, sparsity, order, maxiter, notconv; |
| 64 | + double eps, avg_iter; |
| 65 | + |
| 66 | + void CompareEigen(ModuleBase::ComplexMatrix &psi, double *precondition) |
| 67 | + { |
| 68 | + //calculate eigenvalues by LAPACK; |
| 69 | + double* e_lapack = new double[npw]; |
| 70 | + ModuleBase::ComplexMatrix ev; |
| 71 | + lapackEigen(npw, DIAGOTEST::hmatrix, ev, e_lapack,DETAILINFO); |
| 72 | + |
| 73 | + //do Diago_David::diag() |
| 74 | + double* en = new double[npw]; |
| 75 | + |
| 76 | + Hamilt_PW hpw; |
| 77 | + Diago_David dav(&hpw); |
| 78 | + clock_t start,end; |
| 79 | + start = clock(); |
| 80 | + dav.diag(psi,en,npw,nband,precondition,order,eps,maxiter,notconv,avg_iter); |
| 81 | + end = clock(); |
| 82 | + if (DETAILINFO) std::cout<<"diag Run time: "<<(double)(end - start) / CLOCKS_PER_SEC<<" S, notconv=" |
| 83 | + << notconv << ", avg_iter=" << avg_iter << std::endl; |
| 84 | + |
| 85 | + for(int i=0;i<nband;i++) |
| 86 | + { |
| 87 | + EXPECT_NEAR(en[i],e_lapack[i],CONVTHRESHOLD); |
| 88 | + } |
| 89 | + |
| 90 | + delete [] en; |
| 91 | + delete [] e_lapack; |
| 92 | + } |
| 93 | +}; |
| 94 | + |
| 95 | +class DiagoDavTest : public ::testing::TestWithParam<DiagoDavPrepare> {}; |
| 96 | + |
| 97 | +TEST_P(DiagoDavTest,RandomHamilt) |
| 98 | +{ |
| 99 | + DiagoDavPrepare ddp = GetParam(); |
| 100 | + if (DETAILINFO) std::cout << "npw=" << ddp.npw << ", nband=" << ddp.nband << ", sparsity=" |
| 101 | + << ddp.sparsity << ", eps=" << ddp.eps << std::endl; |
| 102 | + |
| 103 | + HPsi hpsi(ddp.nband,ddp.npw,ddp.sparsity); |
| 104 | + DIAGOTEST::hmatrix = hpsi.hamilt(); |
| 105 | + DIAGOTEST::npw = ddp.npw; |
| 106 | + ModuleBase::ComplexMatrix psi = hpsi.psi(); |
| 107 | + |
| 108 | + ddp.CompareEigen(psi,hpsi.precond()); |
| 109 | +} |
| 110 | + |
| 111 | + |
| 112 | +INSTANTIATE_TEST_SUITE_P(VerifyDiag,DiagoDavTest,::testing::Values( |
| 113 | + //DiagoDavPrepare(int nband, int npw, int sparsity, int order,double eps,int maxiter) |
| 114 | + DiagoDavPrepare(10,100,0,4,1e-5,500), |
| 115 | + DiagoDavPrepare(10,100,5,4,1e-5,500), |
| 116 | + DiagoDavPrepare(10,100,9,4,1e-5,500), |
| 117 | + DiagoDavPrepare(20,500,0,4,1e-5,500), |
| 118 | + DiagoDavPrepare(20,500,5,4,1e-5,500), |
| 119 | + DiagoDavPrepare(20,500,9,4,1e-5,500), |
| 120 | + DiagoDavPrepare(10,1000,8,4,1e-5,500) |
| 121 | + //DiagoDavPrepare(20,2000,8,4,1e-5,500) |
| 122 | +)); |
| 123 | + |
| 124 | +TEST(DiagoDavRealSystemTest,dataH) |
| 125 | +{ |
| 126 | + ModuleBase::ComplexMatrix hmatrix; |
| 127 | + std::ifstream ifs("data-H"); |
| 128 | + DIAGOTEST::readh(ifs,hmatrix); |
| 129 | + DIAGOTEST::hmatrix = hmatrix; |
| 130 | + DIAGOTEST::npw = hmatrix.nc; |
| 131 | + |
| 132 | + DiagoDavPrepare ddp(10,DIAGOTEST::npw,0,2,1e-5,500); |
| 133 | + HPsi hpsi(10,DIAGOTEST::npw); |
| 134 | + ModuleBase::ComplexMatrix psi = hpsi.psi(); |
| 135 | + |
| 136 | + ddp.CompareEigen(psi,hpsi.precond()); |
| 137 | +} |
| 138 | + |
| 139 | +int main(int argc, char **argv) |
| 140 | +{ |
| 141 | + MPI_Init(&argc, &argv); |
| 142 | + testing::InitGoogleTest(&argc, argv); |
| 143 | + int result = RUN_ALL_TESTS(); |
| 144 | + MPI_Finalize(); |
| 145 | + return result; |
| 146 | +} |
0 commit comments