diff --git a/source/module_base/blacs_connector.h b/source/module_base/blacs_connector.h index d29ac47fae..98fd06fcde 100644 --- a/source/module_base/blacs_connector.h +++ b/source/module_base/blacs_connector.h @@ -33,7 +33,8 @@ extern "C" void Cblacs_gridmap(int* icontxt, int *usermap, int ldumap, int nprow, int npcol); // Informational and Miscellaneous void Cblacs_gridinfo(int icontxt, int* nprow, int *npcol, int *myprow, int *mypcol); - int Cblacs_pnum(int icontxt, int prow, int pcol); - void Cblacs_pcoord(int icontxt, int pnum, int *prow, int *pcol); + void Cblacs_gridinit(int* icontxt, char* layout, int nprow, int npcol); + int Cblacs_pnum(int icontxt, int prow, int pcol); + void Cblacs_pcoord(int icontxt, int pnum, int *prow, int *pcol); void Cblacs_exit(int icontxt); } diff --git a/source/module_base/test/CMakeLists.txt b/source/module_base/test/CMakeLists.txt index dc55ca3aed..c281e8beab 100644 --- a/source/module_base/test/CMakeLists.txt +++ b/source/module_base/test/CMakeLists.txt @@ -4,6 +4,11 @@ AddTest( LIBS ${math_libs} SOURCES blas_connector_test.cpp ) +AddTest( + TARGET base_blacs_connector + LIBS ${math_libs} + SOURCES blacs_connector_test.cpp +) AddTest( TARGET base_timer SOURCES timer_test.cpp ../timer.cpp diff --git a/source/module_base/test/blacs_connector_test.cpp b/source/module_base/test/blacs_connector_test.cpp new file mode 100644 index 0000000000..17cd92d268 --- /dev/null +++ b/source/module_base/test/blacs_connector_test.cpp @@ -0,0 +1,53 @@ +#include "../blacs_connector.h" + +#include + +#include "gtest/gtest.h" + +/************************************************ + * unit test of functions in blacs_connector.h + ***********************************************/ + +/** + * - Tested Function + * - Cblacs_gridinit + * - Initializes a grid of processors with a given number of rows and columns. + * The function creates a cartesian topology of all the processors initialized + * by the BLS library. In this topology, each processor is identified by its + * coordinates (row, col) in the grid. + */ + +TEST(blacs_connector, Cblacs_gridinit) +{ + int icontxt; + char layout[] = "ROW"; + int nprow = 1; + int npcol = 1; + + int myid, nprocs; + Cblacs_pinfo(&myid, &nprocs); + Cblacs_get(-1, 0, &icontxt); + + // Call the Cblacs_gridinit() function + Cblacs_gridinit(&icontxt, layout, nprow, npcol); + + // Check if the grid context handle is created successfully + EXPECT_EQ(icontxt, 0); +} + +int main(int argc, char** argv) +{ + int myrank; + int mysize; + + MPI_Init(&argc, &argv); + MPI_Comm_size(MPI_COMM_WORLD, &mysize); + MPI_Comm_rank(MPI_COMM_WORLD, &myrank); + + testing::InitGoogleTest(&argc, argv); + + int result = 0; + result = RUN_ALL_TESTS(); + MPI_Finalize(); + return 0; +} diff --git a/source/module_elecstate/potentials/H_TDDFT_pw.cpp b/source/module_elecstate/potentials/H_TDDFT_pw.cpp index c1b67cc64d..ef47e62b90 100644 --- a/source/module_elecstate/potentials/H_TDDFT_pw.cpp +++ b/source/module_elecstate/potentials/H_TDDFT_pw.cpp @@ -312,19 +312,19 @@ double H_TDDFT_pw::cal_v_time_heaviside() void H_TDDFT_pw::prepare(const UnitCell& cell, int& dir) { - if (dir == 0) + if (dir == 1) { bvec[0] = cell.G.e11; bvec[1] = cell.G.e12; bvec[2] = cell.G.e13; } - else if (dir == 1) + else if (dir == 2) { bvec[0] = cell.G.e21; bvec[1] = cell.G.e22; bvec[2] = cell.G.e23; } - else if (dir == 2) + else if (dir == 3) { bvec[0] = cell.G.e31; bvec[1] = cell.G.e32; diff --git a/source/module_hamilt_lcao/hamilt_lcaodft/FORCE_STRESS.cpp b/source/module_hamilt_lcao/hamilt_lcaodft/FORCE_STRESS.cpp index b72c9e26c1..29251e31da 100644 --- a/source/module_hamilt_lcao/hamilt_lcaodft/FORCE_STRESS.cpp +++ b/source/module_hamilt_lcao/hamilt_lcaodft/FORCE_STRESS.cpp @@ -16,30 +16,32 @@ double Force_Stress_LCAO::force_invalid_threshold_ev = 0.00; double Force_Stress_LCAO::output_acc = 1.0e-8; -Force_Stress_LCAO::Force_Stress_LCAO(Record_adj& ra, const int nat_in) : - RA(&ra), f_pw(nat_in), nat(nat_in){} -Force_Stress_LCAO::~Force_Stress_LCAO() {} - -void Force_Stress_LCAO::getForceStress( - const bool isforce, - const bool isstress, - const bool istestf, - const bool istests, - Local_Orbital_Charge& loc, - const elecstate::ElecState* pelec, - const psi::Psi* psid, - const psi::Psi>* psi, - LCAO_Hamilt &uhm, - ModuleBase::matrix& fcs, - ModuleBase::matrix & scs, - const Structure_Factor& sf, - const K_Vectors& kv, - ModulePW::PW_Basis* rhopw, +Force_Stress_LCAO::Force_Stress_LCAO(Record_adj& ra, const int nat_in) : RA(&ra), f_pw(nat_in), nat(nat_in) +{ +} +Force_Stress_LCAO::~Force_Stress_LCAO() +{ +} + +void Force_Stress_LCAO::getForceStress(const bool isforce, + const bool isstress, + const bool istestf, + const bool istests, + Local_Orbital_Charge& loc, + const elecstate::ElecState* pelec, + const psi::Psi* psid, + const psi::Psi>* psi, + LCAO_Hamilt& uhm, + ModuleBase::matrix& fcs, + ModuleBase::matrix& scs, + const Structure_Factor& sf, + const K_Vectors& kv, + ModulePW::PW_Basis* rhopw, #ifdef __EXX - Exx_LRI& exx_lri_double, - Exx_LRI>& exx_lri_complex, + Exx_LRI& exx_lri_double, + Exx_LRI>& exx_lri_complex, #endif - ModuleSymmetry::Symmetry* symm) + ModuleSymmetry::Symmetry* symm) { ModuleBase::TITLE("Force_Stress_LCAO", "getForceStress"); ModuleBase::timer::tick("Force_Stress_LCAO", "getForceStress"); @@ -154,13 +156,13 @@ void Force_Stress_LCAO::getForceStress( svl_dphi, #endif uhm, - kv); - //implement vdw force or stress here - // Peize Lin add 2014-04-04, update 2021-03-09 - // jiyy add 2019-05-18, update 2021-05-02 - ModuleBase::matrix force_vdw; - ModuleBase::matrix stress_vdw; - + kv); + // implement vdw force or stress here + // Peize Lin add 2014-04-04, update 2021-03-09 + // jiyy add 2019-05-18, update 2021-05-02 + ModuleBase::matrix force_vdw; + ModuleBase::matrix stress_vdw; + auto vdw_solver = vdw::make_vdw(GlobalC::ucell, INPUT); if (vdw_solver != nullptr) { @@ -219,13 +221,13 @@ void Force_Stress_LCAO::getForceStress( if (isforce) { force_dftu.create(nat, 3); - } - if(isstress) - { - stress_dftu.create(3, 3); - } - GlobalC::dftu.force_stress(loc.dm_gamma, loc.dm_k, *uhm.LM, force_dftu, stress_dftu, kv); - } + } + if (isstress) + { + stress_dftu.create(3, 3); + } + GlobalC::dftu.force_stress(loc.dm_gamma, loc.dm_k, *uhm.LM, force_dftu, stress_dftu, kv); + } if (!GlobalV::GAMMA_ONLY_LOCAL) this->flk.finish_k(); @@ -366,35 +368,37 @@ void Force_Stress_LCAO::getForceStress( GlobalC::ld.save_npy_f(fcs, "f_tot.npy", GlobalC::ucell.nat); // Ty/Bohr, F_tot if (GlobalV::deepks_scf) { - GlobalC::ld.save_npy_f(fcs - GlobalC::ld.F_delta, "f_base.npy", GlobalC::ucell.nat); //Ry/Bohr, F_base + GlobalC::ld.save_npy_f(fcs - GlobalC::ld.F_delta, "f_base.npy", GlobalC::ucell.nat); // Ry/Bohr, F_base - if(GlobalV::GAMMA_ONLY_LOCAL) + if (GlobalV::GAMMA_ONLY_LOCAL) { GlobalC::ld.cal_gdmx(loc.dm_gamma[0], - GlobalC::ucell, - GlobalC::ORB, - GlobalC::GridD, - pv->trace_loc_row, - pv->trace_loc_col, - isstress); + GlobalC::ucell, + GlobalC::ORB, + GlobalC::GridD, + pv->trace_loc_row, + pv->trace_loc_col, + isstress); } else - { - GlobalC::ld.cal_gdmx_k(loc.dm_k, - GlobalC::ucell, - GlobalC::ORB, - GlobalC::GridD, - pv->trace_loc_row, - pv->trace_loc_col, - kv.nks, - kv.kvec_d, - isstress); + { + GlobalC::ld.cal_gdmx_k(loc.dm_k, + GlobalC::ucell, + GlobalC::ORB, + GlobalC::GridD, + pv->trace_loc_row, + pv->trace_loc_col, + kv.nks, + kv.kvec_d, + isstress); } - if(GlobalV::deepks_out_unittest) GlobalC::ld.check_gdmx(GlobalC::ucell.nat); + if (GlobalV::deepks_out_unittest) + GlobalC::ld.check_gdmx(GlobalC::ucell.nat); GlobalC::ld.cal_gvx(GlobalC::ucell.nat); - if(GlobalV::deepks_out_unittest) GlobalC::ld.check_gvx(GlobalC::ucell.nat); - GlobalC::ld.save_npy_gvx(GlobalC::ucell.nat);// /Bohr, grad_vx + if (GlobalV::deepks_out_unittest) + GlobalC::ld.check_gvx(GlobalC::ucell.nat); + GlobalC::ld.save_npy_gvx(GlobalC::ucell.nat); // /Bohr, grad_vx } else { @@ -712,7 +716,7 @@ void Force_Stress_LCAO::calForceStressIntegralPart(const bool isGammaOnly, #else ModuleBase::matrix& svl_dphi, #endif - LCAO_Hamilt &uhm, + LCAO_Hamilt& uhm, const K_Vectors& kv) { if (isGammaOnly) @@ -759,7 +763,7 @@ void Force_Stress_LCAO::calForceStressIntegralPart(const bool isGammaOnly, svl_dphi, #endif uhm, - kv); + kv); } return; } diff --git a/source/module_hamilt_lcao/module_tddft/CMakeLists.txt b/source/module_hamilt_lcao/module_tddft/CMakeLists.txt index 8e0adf2b7d..f41dad565d 100644 --- a/source/module_hamilt_lcao/module_tddft/CMakeLists.txt +++ b/source/module_hamilt_lcao/module_tddft/CMakeLists.txt @@ -19,4 +19,8 @@ if(ENABLE_LCAO) add_coverage(tddft) endif() + IF (BUILD_TESTING) + add_subdirectory(test) + endif() + endif() diff --git a/source/module_hamilt_lcao/module_tddft/evolve_psi.cpp b/source/module_hamilt_lcao/module_tddft/evolve_psi.cpp index ebd578efdb..c30f9891c7 100644 --- a/source/module_hamilt_lcao/module_tddft/evolve_psi.cpp +++ b/source/module_hamilt_lcao/module_tddft/evolve_psi.cpp @@ -67,7 +67,7 @@ void evolve_psi(const int nband, /// @input Stmp, Htmp, print_matrix /// @output U_operator Propagator prop(propagator, pv); - prop.compute_propagator(nband, nlocal, Stmp, Htmp, H_laststep, U_operator, print_matrix); + prop.compute_propagator(nlocal, Stmp, Htmp, H_laststep, U_operator, print_matrix); // (3)->>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> diff --git a/source/module_hamilt_lcao/module_tddft/propagator.cpp b/source/module_hamilt_lcao/module_tddft/propagator.cpp index f2711ed2f8..c48e4142ff 100644 --- a/source/module_hamilt_lcao/module_tddft/propagator.cpp +++ b/source/module_hamilt_lcao/module_tddft/propagator.cpp @@ -23,8 +23,7 @@ inline int globalIndex(int localindex, int nblk, int nprocs, int myproc) return gIndex; } -void Propagator::compute_propagator(const int nband, - const int nlocal, +void Propagator::compute_propagator(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, const std::complex* H_laststep, @@ -35,16 +34,16 @@ void Propagator::compute_propagator(const int nband, switch (ptype) { case 0: - compute_propagator_cn2(nband, nlocal, Stmp, Htmp, U_operator, print_matrix); + compute_propagator_cn2(nlocal, Stmp, Htmp, U_operator, print_matrix); break; case 1: tag = 1; - compute_propagator_taylor(nband, nlocal, Stmp, Htmp, U_operator, print_matrix, tag); + compute_propagator_taylor(nlocal, Stmp, Htmp, U_operator, print_matrix, tag); break; case 2: - compute_propagator_etrs(nband, nlocal, Stmp, Htmp, H_laststep, U_operator, print_matrix); + compute_propagator_etrs(nlocal, Stmp, Htmp, H_laststep, U_operator, print_matrix); break; @@ -54,8 +53,7 @@ void Propagator::compute_propagator(const int nband, } } -void Propagator::compute_propagator_cn2(const int nband, - const int nlocal, +void Propagator::compute_propagator_cn2(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, std::complex* U_operator, @@ -259,8 +257,7 @@ void Propagator::compute_propagator_cn2(const int nband, delete[] ipiv; } -void Propagator::compute_propagator_taylor(const int nband, - const int nlocal, +void Propagator::compute_propagator_taylor(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, std::complex* U_operator, @@ -580,8 +577,7 @@ void Propagator::compute_propagator_taylor(const int nband, delete[] ipiv; } -void Propagator::compute_propagator_etrs(const int nband, - const int nlocal, +void Propagator::compute_propagator_etrs(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, const std::complex* H_laststep, @@ -593,8 +589,8 @@ void Propagator::compute_propagator_etrs(const int nband, ModuleBase::GlobalFunc::ZEROS(U1, this->ParaV->nloc); ModuleBase::GlobalFunc::ZEROS(U2, this->ParaV->nloc); int tag = 2; - compute_propagator_taylor(nband, nlocal, Stmp, Htmp, U1, print_matrix, tag); - compute_propagator_taylor(nband, nlocal, Stmp, H_laststep, U2, print_matrix, tag); + compute_propagator_taylor(nlocal, Stmp, Htmp, U1, print_matrix, tag); + compute_propagator_taylor(nlocal, Stmp, H_laststep, U2, print_matrix, tag); ScalapackConnector::gemm('N', 'N', nlocal, diff --git a/source/module_hamilt_lcao/module_tddft/propagator.h b/source/module_hamilt_lcao/module_tddft/propagator.h index 7d9e65b2c4..a156ca1355 100644 --- a/source/module_hamilt_lcao/module_tddft/propagator.h +++ b/source/module_hamilt_lcao/module_tddft/propagator.h @@ -24,7 +24,6 @@ class Propagator /** * @brief compute propagator * - * @param[in] nband number of bands * @param[in] nlocal number of orbitals * @param[in] Stmp overlap matrix * @param[in] Htmp H(t+dt/2) or H(t+dt) @@ -32,8 +31,7 @@ class Propagator * @param[in] print_matirx print internal matrix or not * @param[out] U_operator operator of propagator */ - void compute_propagator(const int nband, - const int nlocal, + void compute_propagator(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, const std::complex* H_laststep, @@ -50,15 +48,13 @@ class Propagator /** * @brief compute propagator of method Crank-Nicolson * - * @param[in] nband number of bands * @param[in] nlocal number of orbitals * @param[in] Stmp overlap matrix * @param[in] Htmp H(t+dt/2) or H(t+dt) * @param[in] print_matirx print internal matrix or not * @param[out] U_operator operator of propagator */ - void compute_propagator_cn2(const int nband, - const int nlocal, + void compute_propagator_cn2(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, std::complex* U_operator, @@ -67,7 +63,6 @@ class Propagator /** * @brief compute propagator of method 4th Taylor * - * @param[in] nband number of bands * @param[in] nlocal number of orbitals * @param[in] Stmp overlap matrix * @param[in] Htmp H(t+dt/2) or H(t+dt) @@ -75,8 +70,7 @@ class Propagator * @param[in] tag a parametre different for 4th Taylor and ETRS * @param[out] U_operator operator of propagator */ - void compute_propagator_taylor(const int nband, - const int nlocal, + void compute_propagator_taylor(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, std::complex* U_operator, @@ -86,7 +80,6 @@ class Propagator /** * @brief compute propagator of method ETRS * - * @param[in] nband number of bands * @param[in] nlocal number of orbitals * @param[in] Stmp overlap matrix * @param[in] Htmp H(t+dt/2) or H(t+dt) @@ -94,8 +87,7 @@ class Propagator * @param[in] print_matirx print internal matrix or not * @param[out] U_operator operator of propagator */ - void compute_propagator_etrs(const int nband, - const int nlocal, + void compute_propagator_etrs(const int nlocal, const std::complex* Stmp, const std::complex* Htmp, const std::complex* H_laststep, diff --git a/source/module_hamilt_lcao/module_tddft/test/CMakeLists.txt b/source/module_hamilt_lcao/module_tddft/test/CMakeLists.txt new file mode 100644 index 0000000000..675806ca8d --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/CMakeLists.txt @@ -0,0 +1,34 @@ +remove_definitions(-D __MPI) + +add_library(tddft_test_lib tddft_test.cpp) + +AddTest( + TARGET tddft_middle_hamilt_test + LIBS ${math_libs} base device tddft_test_lib + SOURCES middle_hamilt_test.cpp ../middle_hamilt.cpp +) + +AddTest( + TARGET tddft_bandenergy_test + LIBS ${math_libs} base device tddft_test_lib + SOURCES bandenergy_test.cpp ../bandenergy.cpp +) + +AddTest( + TARGET tddft_norm_psi_test + LIBS ${math_libs} base device tddft_test_lib + SOURCES norm_psi_test.cpp ../norm_psi.cpp +) + +AddTest( + TARGET tddft_upsi_test + LIBS ${math_libs} base device tddft_test_lib + SOURCES upsi_test1.cpp upsi_test2.cpp upsi_test3.cpp ../upsi.cpp +) + +AddTest( + TARGET tddft_propagator_test + LIBS ${math_libs} base device tddft_test_lib + SOURCES propagator_test1.cpp propagator_test2.cpp propagator_test3.cpp ../propagator.cpp +) + diff --git a/source/module_hamilt_lcao/module_tddft/test/bandenergy_test.cpp b/source/module_hamilt_lcao/module_tddft/test/bandenergy_test.cpp new file mode 100644 index 0000000000..45ef2d6359 --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/bandenergy_test.cpp @@ -0,0 +1,109 @@ +#include "module_hamilt_lcao/module_tddft/bandenergy.h" + +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "module_hamilt_lcao/module_tddft/evolve_elec.h" +#include "tddft_test.h" + +/************************************************ + * unit test of functions in bandenergy.h + ***********************************************/ + +/** + * - Tested Function + * - compute_ekb + * - compute band energy ekb . + */ + +#define doublethreshold 1e-8 +double module_tddft::Evolve_elec::td_print_eij = -1; + +Parallel_Orbitals::Parallel_Orbitals() +{ +} +Parallel_Orbitals::~Parallel_Orbitals() +{ +} + +TEST(BandEnergyTest, testBandEnergy) +{ + std::complex* psi_k; + std::complex* Htmp; + double* ekb; + int nband = 3; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nlocal; + pv->nloc_wfc = nlocal * nband; + pv->ncol = nlocal; + pv->nrow = nlocal; + pv->ncol_bands = nband; + pv->dim0 = 1; + pv->dim1 = 1; + pv->nb = 1; + + int dim[2]; + int period[2] = {1, 1}; + int reorder = 0; + dim[0] = nprow; + dim[1] = npcol; + + MPI_Cart_create(MPI_COMM_WORLD, 2, dim, period, reorder, &pv->comm_2D); + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nband, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow), + lld1 = numroc_(&nband, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + descinit_(pv->desc_wfc, &nlocal, &nband, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + descinit_(pv->desc_Eij, &nband, &nband, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + Htmp = new std::complex[nlocal * nlocal]; + psi_k = new std::complex[nlocal * nband]; + ekb = new double[nband]; + + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nlocal; ++j) + { + if (i == j) + { + Htmp[i * nlocal + j] = std::complex(1.0, 0.0); + } + } + } + Htmp[1] = 0.5; + Htmp[4] = 0.5; + + psi_k[0] = 1.0; + psi_k[1] = 1.0; + psi_k[2] = 0.0; + psi_k[3] = 0.0; + psi_k[4] = 2.0; + psi_k[5] = 1.0; + psi_k[6] = 1.0; + psi_k[7] = 0.0; + psi_k[8] = 3.0; + psi_k[9] = 0.0; + psi_k[10] = 0.0; + psi_k[11] = 1.0; + + // Call the function + module_tddft::compute_ekb(pv, nband, nlocal, Htmp, psi_k, ekb); + + // Check the results + EXPECT_NEAR(ekb[0], 3.0, doublethreshold); + EXPECT_NEAR(ekb[1], 8.0, doublethreshold); + EXPECT_NEAR(ekb[2], 10.0, doublethreshold); + + delete[] psi_k; + delete[] Htmp; + delete[] ekb; +} diff --git a/source/module_hamilt_lcao/module_tddft/test/middle_hamilt_test.cpp b/source/module_hamilt_lcao/module_tddft/test/middle_hamilt_test.cpp new file mode 100644 index 0000000000..8387a82de7 --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/middle_hamilt_test.cpp @@ -0,0 +1,107 @@ +#include "module_hamilt_lcao/module_tddft/middle_hamilt.h" + +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "tddft_test.h" + +/************************************************ + * unit test of functions in middle_hamilt.h + ***********************************************/ + +/** + * - Tested Function + * - half_Hmatrix + * - compute H(t+dt/2). + */ + +#define doublethreshold 1e-8 + +Parallel_Orbitals::Parallel_Orbitals() +{ +} +Parallel_Orbitals::~Parallel_Orbitals() +{ +} + +TEST(MiddleHamiltTest, testMiddleHamilt) +{ + std::complex* Htmp; + std::complex* Hlaststep; + int nband = 3; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nlocal; + pv->ncol = nlocal; + pv->nrow = nlocal; + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nlocal, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + Htmp = new std::complex[nlocal * nlocal]; + Hlaststep = new std::complex[nlocal * nlocal]; + + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nlocal; ++j) + { + if (i == j) + { + Htmp[i * nlocal + j] = std::complex(1.0, 0.0); + Hlaststep[i * nlocal + j] = std::complex(1.0 + 0.2 * (i * nlocal + j), 0.0); + } + else + { + Hlaststep[i * nlocal + j] = std::complex(0.2 * (i * nlocal + j), 0.0); + } + } + } + + // Call the function + module_tddft::half_Hmatrix(pv, nband, nlocal, Htmp, Hlaststep, print_matrix); + + // Check the results + EXPECT_NEAR(Htmp[0].real(), 1.0, doublethreshold); + EXPECT_NEAR(Htmp[0].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[1].real(), 0.1, doublethreshold); + EXPECT_NEAR(Htmp[1].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[2].real(), 0.2, doublethreshold); + EXPECT_NEAR(Htmp[2].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[3].real(), 0.3, doublethreshold); + EXPECT_NEAR(Htmp[3].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[4].real(), 0.4, doublethreshold); + EXPECT_NEAR(Htmp[4].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[5].real(), 1.5, doublethreshold); + EXPECT_NEAR(Htmp[5].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[6].real(), 0.6, doublethreshold); + EXPECT_NEAR(Htmp[6].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[7].real(), 0.7, doublethreshold); + EXPECT_NEAR(Htmp[7].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[8].real(), 0.8, doublethreshold); + EXPECT_NEAR(Htmp[8].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[9].real(), 0.9, doublethreshold); + EXPECT_NEAR(Htmp[9].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[10].real(), 2.0, doublethreshold); + EXPECT_NEAR(Htmp[10].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[11].real(), 1.1, doublethreshold); + EXPECT_NEAR(Htmp[11].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[12].real(), 1.2, doublethreshold); + EXPECT_NEAR(Htmp[12].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[13].real(), 1.3, doublethreshold); + EXPECT_NEAR(Htmp[13].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[14].real(), 1.4, doublethreshold); + EXPECT_NEAR(Htmp[14].imag(), 0.0, doublethreshold); + EXPECT_NEAR(Htmp[15].real(), 2.5, doublethreshold); + EXPECT_NEAR(Htmp[15].imag(), 0.0, doublethreshold); + + delete[] Htmp; + delete[] Hlaststep; +} diff --git a/source/module_hamilt_lcao/module_tddft/test/norm_psi_test.cpp b/source/module_hamilt_lcao/module_tddft/test/norm_psi_test.cpp new file mode 100644 index 0000000000..5a1bdc0361 --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/norm_psi_test.cpp @@ -0,0 +1,125 @@ +#include "module_hamilt_lcao/module_tddft/norm_psi.h" + +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "tddft_test.h" + +/************************************************ + * unit test of functions in norm_psi.h + ***********************************************/ + +/** + * - Tested Function + * - norm_psi + * - normalize the wave function. + */ + +#define doublethreshold 1e-8 + +Parallel_Orbitals::Parallel_Orbitals() +{ +} +Parallel_Orbitals::~Parallel_Orbitals() +{ +} + +TEST(NormPsiTest, testNormPsi) +{ + std::complex* psi_k; + std::complex* Stmp; + int nband = 3; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nlocal; + pv->nloc_wfc = nlocal * nband; + pv->ncol = nlocal; + pv->nrow = nlocal; + pv->ncol_bands = nband; + pv->dim0 = 1; + pv->dim1 = 1; + pv->nb = 1; + + int dim[2]; + int period[2] = {1, 1}; + int reorder = 0; + dim[0] = nprow; + dim[1] = npcol; + + MPI_Cart_create(MPI_COMM_WORLD, 2, dim, period, reorder, &pv->comm_2D); + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nband, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow), + lld1 = numroc_(&nband, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + descinit_(pv->desc_wfc, &nlocal, &nband, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + descinit_(pv->desc_Eij, &nband, &nband, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + Stmp = new std::complex[nlocal * nlocal]; + psi_k = new std::complex[nlocal * nband]; + + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nlocal; ++j) + { + if (i == j) + { + Stmp[i * nlocal + j] = std::complex(1.0, 0.0); + } + } + } + Stmp[1] = 0.5; + Stmp[4] = 0.5; + + psi_k[0] = 1.0; + psi_k[1] = 1.0; + psi_k[2] = 0.0; + psi_k[3] = 0.0; + psi_k[4] = 2.0; + psi_k[5] = 1.0; + psi_k[6] = 1.0; + psi_k[7] = 0.0; + psi_k[8] = 3.0; + psi_k[9] = 0.0; + psi_k[10] = 0.0; + psi_k[11] = 1.0; + + // Call the function + module_tddft::norm_psi(pv, nband, nlocal, Stmp, psi_k, print_matrix); + + // Check the results + EXPECT_NEAR(psi_k[0].real(), 0.577350269189626, doublethreshold); + EXPECT_NEAR(psi_k[0].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[1].real(), 0.577350269189626, doublethreshold); + EXPECT_NEAR(psi_k[1].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[2].real(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[2].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[3].real(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[3].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[4].real(), 0.707106781186547, doublethreshold); + EXPECT_NEAR(psi_k[4].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[5].real(), 0.353553390593274, doublethreshold); + EXPECT_NEAR(psi_k[5].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[6].real(), 0.353553390593274, doublethreshold); + EXPECT_NEAR(psi_k[6].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[7].real(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[7].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[8].real(), 0.948683298050514, doublethreshold); + EXPECT_NEAR(psi_k[8].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[9].real(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[9].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[10].real(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[10].imag(), 0.0, doublethreshold); + EXPECT_NEAR(psi_k[11].real(), 0.316227766016838, doublethreshold); + EXPECT_NEAR(psi_k[11].imag(), 0.0, doublethreshold); + + delete[] psi_k; + delete[] Stmp; +} diff --git a/source/module_hamilt_lcao/module_tddft/test/propagator_test1.cpp b/source/module_hamilt_lcao/module_tddft/test/propagator_test1.cpp new file mode 100644 index 0000000000..ce17a13185 --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/propagator_test1.cpp @@ -0,0 +1,114 @@ +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "module_hamilt_lcao/module_tddft/propagator.h" +#include "module_io/input.h" +#include "tddft_test.h" + +/************************************************ + * unit test of functions in propagator.h + ***********************************************/ + +/** + * - Tested Function + * - Propagator::compute_propagator_cn2 + * - compute propagator of method Crank-Nicolson. + */ + +Input INPUT; +#define doublethreshold 1e-8 + +Parallel_Orbitals::Parallel_Orbitals() +{ +} +Parallel_Orbitals::~Parallel_Orbitals() +{ +} + +TEST(PropagatorTest, testPropagatorCN) +{ + std::complex* U_operator; + std::complex* Stmp; + std::complex* Htmp; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nlocal; + pv->ncol = nlocal; + INPUT.mdp.md_dt = 4; + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nlocal, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + U_operator = new std::complex[nlocal * nlocal]; + Stmp = new std::complex[nlocal * nlocal]; + Htmp = new std::complex[nlocal * nlocal]; + + for (int i = 0; i < nlocal * nlocal; ++i) + { + U_operator[i] = std::complex(0.0, 0.0); + } + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nlocal; ++j) + { + if (i == j) + { + Htmp[i * nlocal + j] = std::complex(1.0, 0.0); + Stmp[i * nlocal + j] = std::complex(1.0, 0.0); + } + } + } + Stmp[1] = 0.5; + Stmp[4] = 0.5; + + // Call the function + int propagator = 0; + module_tddft::Propagator prop(propagator, pv); + prop.compute_propagator(nlocal, Stmp, Htmp, nullptr, U_operator, print_matrix); + + // Check the results + EXPECT_NEAR(U_operator[0].real(), -0.107692307692308, doublethreshold); + EXPECT_NEAR(U_operator[0].imag(), -0.861538461538462, doublethreshold); + EXPECT_NEAR(U_operator[1].real(), 0.492307692307692, doublethreshold); + EXPECT_NEAR(U_operator[1].imag(), -0.0615384615384615, doublethreshold); + EXPECT_NEAR(U_operator[2].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[2].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[3].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[3].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[4].real(), 0.492307692307692, doublethreshold); + EXPECT_NEAR(U_operator[4].imag(), -0.0615384615384615, doublethreshold); + EXPECT_NEAR(U_operator[5].real(), -0.107692307692308, doublethreshold); + EXPECT_NEAR(U_operator[5].imag(), -0.861538461538462, doublethreshold); + EXPECT_NEAR(U_operator[6].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[6].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[7].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[7].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[8].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[8].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[9].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[9].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[10].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[10].imag(), -1.0, doublethreshold); + EXPECT_NEAR(U_operator[11].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[11].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[12].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[12].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[13].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[13].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[14].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[14].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[15].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[15].imag(), -1.0, doublethreshold); + + delete[] U_operator; + delete[] Htmp; + delete[] Stmp; +} \ No newline at end of file diff --git a/source/module_hamilt_lcao/module_tddft/test/propagator_test2.cpp b/source/module_hamilt_lcao/module_tddft/test/propagator_test2.cpp new file mode 100644 index 0000000000..e015e17f45 --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/propagator_test2.cpp @@ -0,0 +1,119 @@ +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "module_hamilt_lcao/module_tddft/propagator.h" +#include "module_io/input.h" +#include "tddft_test.h" + +/************************************************ + * unit test of functions in propagator.h + ***********************************************/ + +/** + * - Tested Function + * - Propagator::compute_propagator_taylor + * - compute propagator of method 4th Taylor expansion. + */ + +#define doublethreshold 1e-8 + +TEST(PropagatorTest, testPropagatorTaylor) +{ + std::complex* U_operator; + std::complex* Stmp; + std::complex* Htmp; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nlocal; + pv->ncol = nlocal; + pv->nrow = nlocal; + pv->dim0 = 1; + pv->dim1 = 1; + pv->nb = 1; + + int dim[2]; + int period[2] = {1, 1}; + int reorder = 0; + dim[0] = nprow; + dim[1] = npcol; + + MPI_Cart_create(MPI_COMM_WORLD, 2, dim, period, reorder, &pv->comm_2D); + + INPUT.mdp.md_dt = 4; + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nlocal, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + U_operator = new std::complex[nlocal * nlocal]; + Stmp = new std::complex[nlocal * nlocal]; + Htmp = new std::complex[nlocal * nlocal]; + + for (int i = 0; i < nlocal * nlocal; ++i) + { + U_operator[i] = std::complex(0.0, 0.0); + } + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nlocal; ++j) + { + if (i == j) + { + Htmp[i * nlocal + j] = std::complex(1.0, 0.0); + Stmp[i * nlocal + j] = std::complex(1.0, 0.0); + } + } + } + Stmp[1] = 0.5; + Stmp[4] = 0.5; + + // Call the function + int propagator = 1; + module_tddft::Propagator prop(propagator, pv); + prop.compute_propagator(nlocal, Stmp, Htmp, nullptr, U_operator, print_matrix); + + // Check the results + EXPECT_NEAR(U_operator[0].real(), 1.95473251028807, doublethreshold); + EXPECT_NEAR(U_operator[0].imag(), 2.8641975308642, doublethreshold); + EXPECT_NEAR(U_operator[1].real(), -1.7119341563786, doublethreshold); + EXPECT_NEAR(U_operator[1].imag(), -3.80246913580247, doublethreshold); + EXPECT_NEAR(U_operator[2].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[2].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[3].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[3].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[4].real(), -1.7119341563786, doublethreshold); + EXPECT_NEAR(U_operator[4].imag(), -3.80246913580247, doublethreshold); + EXPECT_NEAR(U_operator[5].real(), 1.95473251028807, doublethreshold); + EXPECT_NEAR(U_operator[5].imag(), 2.8641975308642, doublethreshold); + EXPECT_NEAR(U_operator[6].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[6].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[7].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[7].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[8].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[8].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[9].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[9].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[10].real(), -0.333333333333333, doublethreshold); + EXPECT_NEAR(U_operator[10].imag(), -0.666666666666667, doublethreshold); + EXPECT_NEAR(U_operator[11].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[11].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[12].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[12].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[13].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[13].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[14].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[14].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[15].real(), -0.333333333333333, doublethreshold); + EXPECT_NEAR(U_operator[15].imag(), -0.666666666666667, doublethreshold); + + delete[] U_operator; + delete[] Htmp; + delete[] Stmp; +} \ No newline at end of file diff --git a/source/module_hamilt_lcao/module_tddft/test/propagator_test3.cpp b/source/module_hamilt_lcao/module_tddft/test/propagator_test3.cpp new file mode 100644 index 0000000000..4023b518fe --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/propagator_test3.cpp @@ -0,0 +1,124 @@ +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "module_hamilt_lcao/module_tddft/propagator.h" +#include "module_io/input.h" +#include "tddft_test.h" + +/************************************************ + * unit test of functions in propagator.h + ***********************************************/ + +/** + * - Tested Function + * - Propagator::compute_propagator_etrs + * - compute propagator of method ETRS. + */ + +#define doublethreshold 1e-8 + +TEST(PropagatorTest, testPropagatorETRS) +{ + std::complex* U_operator; + std::complex* Stmp; + std::complex* Htmp; + std::complex* Hlaststep; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nlocal; + pv->ncol = nlocal; + pv->nrow = nlocal; + pv->dim0 = 1; + pv->dim1 = 1; + pv->nb = 1; + + int dim[2]; + int period[2] = {1, 1}; + int reorder = 0; + dim[0] = nprow; + dim[1] = npcol; + + MPI_Cart_create(MPI_COMM_WORLD, 2, dim, period, reorder, &pv->comm_2D); + + INPUT.mdp.md_dt = 4; + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nlocal, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + U_operator = new std::complex[nlocal * nlocal]; + Stmp = new std::complex[nlocal * nlocal]; + Htmp = new std::complex[nlocal * nlocal]; + Hlaststep = new std::complex[nlocal * nlocal]; + + for (int i = 0; i < nlocal * nlocal; ++i) + { + U_operator[i] = std::complex(0.0, 0.0); + } + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nlocal; ++j) + { + if (i == j) + { + Htmp[i * nlocal + j] = std::complex(1.0, 0.0); + Hlaststep[i * nlocal + j] = std::complex(1.0, 0.0); + Stmp[i * nlocal + j] = std::complex(1.0, 0.0); + } + } + } + Stmp[1] = 0.5; + Stmp[4] = 0.5; + Hlaststep[0] = 1.1; + + // Call the function + int propagator = 2; + module_tddft::Propagator prop(propagator, pv); + prop.compute_propagator(nlocal, Stmp, Htmp, Hlaststep, U_operator, print_matrix); + + // Check the results + EXPECT_NEAR(U_operator[0].real(), -0.0105865569272976, doublethreshold); + EXPECT_NEAR(U_operator[0].imag(), -0.228336412132297, doublethreshold); + EXPECT_NEAR(U_operator[1].real(), 0.195527023319616, doublethreshold); + EXPECT_NEAR(U_operator[1].imag(), -0.728701844231062, doublethreshold); + EXPECT_NEAR(U_operator[2].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[2].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[3].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[3].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[4].real(), 0.247662246608749, doublethreshold); + EXPECT_NEAR(U_operator[4].imag(), -0.694263679317177, doublethreshold); + EXPECT_NEAR(U_operator[5].real(), -0.0219180003048316, doublethreshold); + EXPECT_NEAR(U_operator[5].imag(), -0.300090839811004, doublethreshold); + EXPECT_NEAR(U_operator[6].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[6].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[7].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[7].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[8].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[8].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[9].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[9].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[10].real(), -0.401041666666667, doublethreshold); + EXPECT_NEAR(U_operator[10].imag(), -0.902777777777778, doublethreshold); + EXPECT_NEAR(U_operator[11].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[11].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[12].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[12].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[13].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[13].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[14].real(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[14].imag(), 0.0, doublethreshold); + EXPECT_NEAR(U_operator[15].real(), -0.401041666666667, doublethreshold); + EXPECT_NEAR(U_operator[15].imag(), -0.902777777777778, doublethreshold); + + delete[] U_operator; + delete[] Htmp; + delete[] Stmp; + delete[] Hlaststep; +} \ No newline at end of file diff --git a/source/module_hamilt_lcao/module_tddft/test/tddft_test.cpp b/source/module_hamilt_lcao/module_tddft/test/tddft_test.cpp new file mode 100644 index 0000000000..fedb46a976 --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/tddft_test.cpp @@ -0,0 +1,50 @@ +#include "tddft_test.h" + +#include +#include + +#include "module_base/blacs_connector.h" +#include "module_base/scalapack_connector.h" +#include "module_basis/module_ao/parallel_orbitals.h" + +/************************************************ + * unit test of module_tddft + ***********************************************/ + +int myprow, nprow, ictxt, mypcol, npcol; + +void MPIInit() +{ + int myrank; + int mysize; + + MPI_Init(NULL, NULL); + MPI_Comm_size(MPI_COMM_WORLD, &mysize); + MPI_Comm_rank(MPI_COMM_WORLD, &myrank); + + // Set up BLACS context and grid + int nprocs; + nprow = 1; + npcol = 1; + Cblacs_pinfo(&myrank, &mysize); + Cblacs_get(-1, 0, &ictxt); + Cblacs_gridinit(&ictxt, "Row", nprow, npcol); + Cblacs_gridinfo(ictxt, &nprow, &npcol, &myprow, &mypcol); +} + +/************************************************ + * unit test of module_tddft + ***********************************************/ +int main(int argc, char** argv) +{ + MPIInit(); + + int result = 0; + testing::InitGoogleTest(&argc, argv); + result = RUN_ALL_TESTS(); + + Cblacs_exit(ictxt); + + // MPI_Finalize(); + return 0; +} diff --git a/source/module_hamilt_lcao/module_tddft/test/tddft_test.h b/source/module_hamilt_lcao/module_tddft/test/tddft_test.h new file mode 100644 index 0000000000..c13f0b0e4e --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/tddft_test.h @@ -0,0 +1,50 @@ +#ifndef __TDDFTTEST +#define __TDDFTTEST +#include + +#include "gtest/gtest.h" +#include "module_basis/module_ao/parallel_orbitals.h" + +using namespace std; +extern int myprow, nprow, ictxt, mypcol, npcol; + +class TDDFTTEST : public testing::Test +{ + public: + static void SetUpTestCase() + { + cout << "\033[32m" + << "============================" + << "\033[0m" << endl; + cout << "\033[32m" + << "= TDDFT MODULE TEST START =" + << "\033[0m" << endl; + cout << "\033[32m" + << "============================" + << "\033[0m" << endl; + } + static void TearDownTestCase() + { + cout << "\033[32m" + << "============================" + << "\033[0m" << endl; + cout << "\033[32m" + << "= TDDFT MODULE TEST END =" + << "\033[0m" << endl; + cout << "\033[32m" + << "============================" + << "\033[0m" << endl; + } + void SetUp() + { + cout << "\033[32m" + << "[ CASE ]" + << "\033[0m" + << " "; + } + void TearDown() + { + } +}; + +#endif diff --git a/source/module_hamilt_lcao/module_tddft/test/upsi_test1.cpp b/source/module_hamilt_lcao/module_tddft/test/upsi_test1.cpp new file mode 100644 index 0000000000..e73221f80b --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/upsi_test1.cpp @@ -0,0 +1,76 @@ +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "module_hamilt_lcao/module_tddft/upsi.h" +#include "tddft_test.h" + +#define doublethreshold 1e-8 + +/************************************************ + * unit test of functions in upsi.h + ***********************************************/ + +/** + * - Tested Function + * - upsi + * - apply U_operator to the wave function of the previous step for new wave function. + */ + +Parallel_Orbitals::Parallel_Orbitals() +{ +} +Parallel_Orbitals::~Parallel_Orbitals() +{ +} + +TEST(UpsiTest, testUpsi1) +{ + std::complex* U_operator; + std::complex* psi_k_laststep; + std::complex* psi_k; + int nband = 2; + int nlocal = 2; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nband; + pv->ncol = nlocal; + pv->ncol_bands = nband; + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nband, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + descinit_(pv->desc_wfc, &nlocal, &nband, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + U_operator = new std::complex[nlocal * nlocal]; + psi_k_laststep = new std::complex[nlocal * nband]; + psi_k = new std::complex[nlocal * nband]; + + for (int i = 0; i < nlocal * nlocal; ++i) + { + U_operator[i] = std::complex(1.0, 0.0); + } + for (int i = 0; i < nlocal * nband; ++i) + { + psi_k_laststep[i] = std::complex(1.0, 0.0); + psi_k[i] = std::complex(0.0, 0.0); + } + + // Call the function + module_tddft::upsi(pv, nband, nlocal, U_operator, psi_k_laststep, psi_k, print_matrix); + + // Check the results + for (int i = 0; i < nlocal * nband; ++i) + { + EXPECT_NEAR(psi_k[i].real(), nlocal, doublethreshold); + EXPECT_NEAR(psi_k[i].imag(), 0.0, doublethreshold); + } + delete[] U_operator; + delete[] psi_k; + delete[] psi_k_laststep; +} diff --git a/source/module_hamilt_lcao/module_tddft/test/upsi_test2.cpp b/source/module_hamilt_lcao/module_tddft/test/upsi_test2.cpp new file mode 100644 index 0000000000..33ca60f7a6 --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/upsi_test2.cpp @@ -0,0 +1,69 @@ +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "module_hamilt_lcao/module_tddft/upsi.h" +#include "tddft_test.h" + +#define doublethreshold 1e-8 + +/************************************************ + * unit test of functions in upsi.h + ***********************************************/ + +/** + * - Tested Function + * - upsi + * - apply U_operator to the wave function of the previous step for new wave function. + */ + +TEST(UpsiTest, testUpsi2) +{ + std::complex* U_operator; + std::complex* psi_k_laststep; + std::complex* psi_k; + int nband = 3; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nband; + pv->ncol = nlocal; + pv->ncol_bands = nband; + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nband, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + descinit_(pv->desc_wfc, &nlocal, &nband, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + U_operator = new std::complex[nlocal * nlocal]; + psi_k_laststep = new std::complex[nlocal * nband]; + psi_k = new std::complex[nlocal * nband]; + + for (int i = 0; i < nlocal * nlocal; ++i) + { + U_operator[i] = std::complex(1.0, 0.0); + } + for (int i = 0; i < nlocal * nband; ++i) + { + psi_k_laststep[i] = std::complex(2.0, 0.0); + psi_k[i] = std::complex(0.0, 0.0); + } + + // Call the function + module_tddft::upsi(pv, nband, nlocal, U_operator, psi_k_laststep, psi_k, print_matrix); + + // Check the results + for (int i = 0; i < nlocal * nband; ++i) + { + EXPECT_NEAR(psi_k[i].real(), 2.0 * nlocal, doublethreshold); + EXPECT_NEAR(psi_k[i].imag(), 0.0, doublethreshold); + } + delete[] U_operator; + delete[] psi_k; + delete[] psi_k_laststep; +} diff --git a/source/module_hamilt_lcao/module_tddft/test/upsi_test3.cpp b/source/module_hamilt_lcao/module_tddft/test/upsi_test3.cpp new file mode 100644 index 0000000000..576adf60ee --- /dev/null +++ b/source/module_hamilt_lcao/module_tddft/test/upsi_test3.cpp @@ -0,0 +1,86 @@ +#include +#include +#include + +#include "module_basis/module_ao/parallel_orbitals.h" +#include "module_hamilt_lcao/module_tddft/upsi.h" +#include "tddft_test.h" + +#define doublethreshold 1e-8 + +/************************************************ + * unit test of functions in upsi.h + ***********************************************/ + +/** + * - Tested Function + * - upsi + * - apply U_operator to the wave function of the previous step for new wave function. + */ + +TEST(UpsiTest, testUpsi3) +{ + std::complex* U_operator; + std::complex* psi_k_laststep; + std::complex* psi_k; + int nband = 3; + int nlocal = 4; + bool print_matrix = false; + Parallel_Orbitals* pv; + pv = new Parallel_Orbitals(); + pv->nloc = nlocal * nband; + pv->ncol = nlocal; + pv->ncol_bands = nband; + + // Initialize input matrices + int info; + int mb = 1, nb = 1, lda = nband, ldc = nlocal; + int irsrc = 0, icsrc = 0, lld = numroc_(&nlocal, &mb, &myprow, &irsrc, &nprow); + descinit_(pv->desc, &nlocal, &nlocal, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + descinit_(pv->desc_wfc, &nlocal, &nband, &mb, &nb, &irsrc, &icsrc, &ictxt, &lld, &info); + + // Initialize data + U_operator = new std::complex[nlocal * nlocal]; + psi_k_laststep = new std::complex[nlocal * nband]; + psi_k = new std::complex[nlocal * nband]; + + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nlocal; ++j) + { + if (i == j) + { + U_operator[i * nlocal + j] = std::complex(1.0, 0.0); + } + else + { + U_operator[i * nlocal + j] = std::complex(0.0, 0.0); + } + } + } + + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nband; ++j) + { + psi_k_laststep[i * nband + j] = std::complex(i * nband + j, 0.0); + psi_k[i * nband + j] = std::complex(0.0, 0.0); + } + } + + // Call the function + module_tddft::upsi(pv, nband, nlocal, U_operator, psi_k_laststep, psi_k, print_matrix); + + // Check the results + for (int i = 0; i < nlocal; ++i) + { + for (int j = 0; j < nband; ++j) + { + EXPECT_NEAR(psi_k[i * nband + j].real(), i * nband + j, doublethreshold); + EXPECT_NEAR(psi_k[i].imag(), 0.0, doublethreshold); + } + } + delete[] U_operator; + delete[] psi_k; + delete[] psi_k_laststep; +}