diff --git a/.github/mrchem-codecov.yml b/.github/mrchem-codecov.yml index 611adbfba..2e86c7383 100644 --- a/.github/mrchem-codecov.yml +++ b/.github/mrchem-codecov.yml @@ -10,3 +10,4 @@ dependencies: - ninja - nlohmann_json - xcfun + - libxc diff --git a/.github/mrchem-gha.yml b/.github/mrchem-gha.yml index 38c86330b..9529b5582 100644 --- a/.github/mrchem-gha.yml +++ b/.github/mrchem-gha.yml @@ -9,3 +9,4 @@ dependencies: - ninja - nlohmann_json - xcfun + - libxc diff --git a/cmake/custom/main.cmake b/cmake/custom/main.cmake index 5cc4ee481..6e8d7fe8b 100644 --- a/cmake/custom/main.cmake +++ b/cmake/custom/main.cmake @@ -56,6 +56,7 @@ set(_build_type ${CMAKE_BUILD_TYPE}) # Order IS important here! include(${PROJECT_SOURCE_DIR}/external/upstream/fetch_nlohmann_json.cmake) include(${PROJECT_SOURCE_DIR}/external/upstream/fetch_xcfun.cmake) +include(${PROJECT_SOURCE_DIR}/external/upstream/fetch_libxc.cmake) include(${PROJECT_SOURCE_DIR}/external/upstream/fetch_eigen3.cmake) include(${PROJECT_SOURCE_DIR}/external/upstream/fetch_mrcpp.cmake) # reset CMAKE_BUILD_TYPE to whatever it was for MRChem diff --git a/doc/users/user_inp.rst b/doc/users/user_inp.rst index fe8b59989..d3760494e 100644 --- a/doc/users/user_inp.rst +++ b/doc/users/user_inp.rst @@ -239,7 +239,7 @@ Here we specify the exchange-correlation functional used in DFT DFT { spin = false # Use spin-polarized functionals - density_cutoff = 0.0 # Cutoff to set XC potential to zero + density_cutoff = 1e-11 # Cutoff to set XC potential to zero $functionals 1.0 # Functional name and coefficient 1.0 # Functional name and coefficient diff --git a/external/upstream/fetch_eigen3.cmake b/external/upstream/fetch_eigen3.cmake index 68a6683c8..be8f9e0bc 100644 --- a/external/upstream/fetch_eigen3.cmake +++ b/external/upstream/fetch_eigen3.cmake @@ -1,4 +1,4 @@ -find_package(Eigen3 CONFIG QUIET +find_package(Eigen3 3.4 CONFIG QUIET NO_CMAKE_PATH NO_CMAKE_PACKAGE_REGISTRY ) diff --git a/external/upstream/fetch_libxc.cmake b/external/upstream/fetch_libxc.cmake new file mode 100644 index 000000000..ecc49c9a5 --- /dev/null +++ b/external/upstream/fetch_libxc.cmake @@ -0,0 +1,36 @@ +find_package(Libxc QUIET CONFIG) + +if(TARGET Libxc::xc) + get_target_property(_loc Libxc::xc LOCATION) + message(STATUS "Libxc::xc: ${_loc} (found version ${Libxc_VERSION})") +else() + message(STATUS "Suitable LibXC could not be located. Fetching and building!") + include(FetchContent) + FetchContent_Declare(libxc_sources + GIT_REPOSITORY https://gitlab.com/libxc/libxc + GIT_TAG 7.0.0 + ) + set(CMAKE_INSTALL_INCLUDEDIR "include/" CACHE STRING "" FORCE) # Creates Libxc subdir in install + set(BUILD_TESTING OFF CACHE BOOL "Build LibXC tests" FORCE) + set(ENABLE_TESTS OFF CACHE BOOL "Enable LibXC tests" FORCE) + set(BUILD_SHARED_LIBS ON CACHE BOOL "Build LibXC shared libs" FORCE) + set(ENABLE_FORTRAN OFF CACHE BOOL "Build LibXC Fortran bindings" FORCE) + set(DISABLE_FXC ON CACHE BOOL "Disable 2nd derivatives (Fxc)" FORCE) + set(DISABLE_KXC ON CACHE BOOL "Disable 3rd derivatives (Kxc)" FORCE) + set(DISABLE_LXC ON CACHE BOOL "Disable 4th derivatives (Lxc)" FORCE) + + FetchContent_MakeAvailable(libxc_sources) + + if(TARGET xc) + if(NOT TARGET Libxc::xc) + add_library(Libxc::xc ALIAS xc) + endif() + + target_include_directories(xc + INTERFACE + $ + $ + $ + ) + endif() +endif() \ No newline at end of file diff --git a/python/mrchem/helpers.py b/python/mrchem/helpers.py index 2404db378..5c6967f15 100644 --- a/python/mrchem/helpers.py +++ b/python/mrchem/helpers.py @@ -122,6 +122,10 @@ def write_scf_fock(user_dict, wf_dict, origin): }, } + # Exchange-Correlation library + if wf_dict["method_type"] in ["dft"]: + fock_dict["xc_library"] = user_dict["DFT"]["xc_library"], + # External electric field if len(user_dict["ExternalFields"]["electric_field"]) > 0: fock_dict["external_operator"] = { @@ -249,6 +253,8 @@ def write_scf_guess(user_dict, wf_dict): "file_CUBE_p": f"{vector_dir}CUBE_p_vector.json", "file_CUBE_a": f"{vector_dir}CUBE_a_vector.json", "file_CUBE_b": f"{vector_dir}CUBE_b_vector.json", + "xc_library": user_dict["DFT"]["xc_library"], + "cutoff": user_dict["DFT"]["density_cutoff"], } return guess_dict @@ -279,6 +285,7 @@ def write_scf_solver(user_dict, wf_dict): "energy_thrs": scf_dict["energy_thrs"], "orbital_thrs": scf_dict["orbital_thrs"], "helmholtz_prec": user_dict["Precisions"]["helmholtz_prec"], + "xc_library": user_dict["DFT"]["xc_library"], "deltascf_method": scf_dict["deltascf_method"], } @@ -466,6 +473,12 @@ def write_rsp_fock(user_dict, wf_dict): }, } + # Exchange-Correlation library + if wf_dict["method_type"] in ["dft"]: + fock_dict["xc_library"] = { + "xc_library": user_dict["DFT"]["xc_library"], + } + # Reaction if user_dict["WaveFunction"]["environment"].lower() != "none": fock_dict["reaction_operator"] = _reaction_operator_handler(user_dict, rsp=True) @@ -496,6 +509,8 @@ def write_rsp_solver(user_dict, wf_dict, d): "property_thrs": user_dict["Response"]["property_thrs"], "helmholtz_prec": user_dict["Precisions"]["helmholtz_prec"], "orth_prec": 1.0e-14, + "xc_library": user_dict["DFT"]["xc_library"], + "cutoff": user_dict["DFT"]["density_cutoff"], } return solver_dict diff --git a/python/mrchem/input_parser/api.py b/python/mrchem/input_parser/api.py index 9ae98a521..cc4831ff7 100644 --- a/python/mrchem/input_parser/api.py +++ b/python/mrchem/input_parser/api.py @@ -338,7 +338,7 @@ def stencil() -> JSONDict: 'name': 'azora_potential_path', 'type': 'str'}], 'name': 'ZORA'}, - { 'keywords': [ { 'default': 0.0, + { 'keywords': [ { 'default': 1e-11, 'name': 'density_cutoff', 'type': 'float'}, { 'default': ' ', @@ -346,7 +346,14 @@ def stencil() -> JSONDict: 'type': 'str'}, { 'default': "not(user['WaveFunction']['restricted'])", 'name': 'spin', - 'type': 'bool'}], + 'type': 'bool'}, + # { 'default': True, + # 'name': 'libxc', + # 'type': 'bool'} + { 'default': 'xcfun', + 'name': 'xc_library', + 'type': 'str'} + ], 'name': 'DFT'}, { 'keywords': [ { 'default': True, 'name': 'dipole_moment', diff --git a/python/template.yml b/python/template.yml index 082bc1327..2a6dce2e9 100644 --- a/python/template.yml +++ b/python/template.yml @@ -368,6 +368,11 @@ sections: default: 0.0 docstring: | Hard cutoff for passing density values to XCFun. + - name: xc_library + type: str + default: xcfun + docstring: | + Runs Libxc or XCFun - name: functionals type: str default: ' ' diff --git a/setup b/setup index 01abb222b..c790a033f 100755 --- a/setup +++ b/setup @@ -47,6 +47,7 @@ def gen_cmake_command(options, arguments): command.append('-DENABLE_MPI={0}'.format(arguments['--mpi'])) command.append('-DENABLE_OPENMP={0}'.format(arguments['--omp'])) command.append('-DCMAKE_BUILD_TYPE={0}'.format(arguments['--type'])) + command.append('-DCMAKE_POLICY_VERSION_MINIMUM=3.5') command.append('-G"{0}"'.format(arguments['--generator'])) if arguments['--cmake-options'] != "''": command.append(arguments['--cmake-options']) diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index aec3eae2b..e2e2bd3fd 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -51,6 +51,7 @@ target_link_libraries(mrchem PRIVATE Eigen3::Eigen PUBLIC + Libxc::xc XCFun::xcfun MRCPP::mrcpp nlohmann_json::nlohmann_json diff --git a/src/driver.cpp b/src/driver.cpp index 454f6346a..42dc50a0c 100644 --- a/src/driver.cpp +++ b/src/driver.cpp @@ -26,6 +26,7 @@ #include #include #include +#include // libxc debug #include "driver.h" #include @@ -248,6 +249,15 @@ json driver::scf::run(const json &json_scf, Molecule &mol) { json json_out = {{"success", true}}; if (json_scf.contains("properties")) driver::init_properties(json_scf["properties"], mol); + /////////////////////////////////////////////////////////// + ////////////////// Setting XC Library ////////////////// + /////////////////////////////////////////////////////////// + std::string xc_lib; + + if (json_scf["fock_operator"].contains("xc_library")) { + xc_lib = json_scf["fock_operator"]["xc_library"][0].get(); + } else {xc_lib = "xcfun";} + /////////////////////////////////////////////////////////// //////////////// Building Fock Operator /////////////// /////////////////////////////////////////////////////////// @@ -312,6 +322,7 @@ json driver::scf::run(const json &json_scf, Molecule &mol) { solver.setRotation(rotation); solver.setLocalize(localize); solver.setMethodName(method); + solver.setLibxc((xc_lib == "libxc") ? true : false); solver.setRelativityName(relativity); solver.setEnvironmentName(environment); solver.setExternalFieldName(external_field); @@ -507,14 +518,18 @@ bool driver::scf::guess_energy(const json &json_guess, Molecule &mol, FockBuilde auto external_field = json_guess["external_field"]; auto localize = json_guess["localize"]; auto rotate = json_guess["rotate"]; + std::string xc_lib = json_guess["xc_library"].get(); + auto cutoff = json_guess["cutoff"]; mrcpp::print::separator(0, '~'); print_utils::text(0, "Calculation ", "Compute initial energy"); print_utils::text(0, "Method ", method); + print_utils::text(0, "XC Library ", (xc_lib == "libxc") ? "LibXC" : "XCFun"); print_utils::text(0, "Relativity ", relativity); print_utils::text(0, "Environment ", environment); print_utils::text(0, "External fields", external_field); print_utils::text(0, "Precision ", print_utils::dbl_to_str(prec, 5, true)); + print_utils::text(0, "Density cutoff ", print_utils::dbl_to_str(cutoff, 5, true)); print_utils::text(0, "Localization ", (localize) ? "On" : "Off"); mrcpp::print::separator(0, '~', 2); @@ -1343,9 +1358,21 @@ void driver::build_fock_operator(const json &json_fock, Molecule &mol, FockBuild auto xc_cutoff = json_xcfunc["cutoff"]; auto xc_funcs = json_xcfunc["functionals"]; auto xc_order = order + 1; + // TODO: Look over and input parser so this is not necessary + std::string xc_lib; + if (json_fock.contains("xc_library")) { + if(json_fock["xc_library"].is_array()){ + xc_lib = json_fock["xc_library"][0].get(); + }else{ + xc_lib = json_fock["xc_library"]["xc_library"].get(); + } + }else{ + xc_lib = "xcfun"; + } mrdft::Factory xc_factory(*MRA); xc_factory.setSpin(xc_spin); + xc_factory.setLibxc((xc_lib == "libxc") ? true : false); xc_factory.setOrder(xc_order); xc_factory.setDensityCutoff(xc_cutoff); for (const auto &f : xc_funcs) { @@ -1354,6 +1381,9 @@ void driver::build_fock_operator(const json &json_fock, Molecule &mol, FockBuild xc_factory.setFunctional(name, coef); } auto mrdft_p = xc_factory.build(); + + mrdft_p->functional().print_functional_references(); + exx = mrdft_p->functional().amountEXX(); if (order == 0) { diff --git a/src/initial_guess/sad.cpp b/src/initial_guess/sad.cpp index 94fc0c66a..2bdf53dad 100644 --- a/src/initial_guess/sad.cpp +++ b/src/initial_guess/sad.cpp @@ -72,6 +72,7 @@ bool initial_guess::sad::setup(OrbitalVector &Phi, double prec, double screen, c print_utils::text(0, "Screening ", print_utils::dbl_to_str(screen, 5, true) + " StdDev"); print_utils::text(0, "Restricted ", (restricted) ? "True" : "False"); print_utils::text(0, "Functional ", "LDA (SVWN5)"); + print_utils::text(0, "XC Library ", (mrdft::Factory::libxc) ? "LibXC" : "XCFun"); print_utils::text(0, "AO basis ", "Hydrogenic orbitals"); print_utils::text(0, "Zeta quality", std::to_string(zeta)); mrcpp::print::separator(0, '~', 2); @@ -90,8 +91,8 @@ bool initial_guess::sad::setup(OrbitalVector &Phi, double prec, double screen, c MomentumOperator p(D_p); NuclearOperator V_nuc(nucs, prec); CoulombOperator J(P_p); - XCOperator XC(mrdft_p); - RankZeroOperator V = V_nuc + J + XC; + XCOperator XC_(mrdft_p); + RankZeroOperator V = V_nuc + J + XC_; auto plevel = Printer::getPrintLevel(); if (plevel == 1) mrcpp::print::header(1, "SAD Initial Guess"); @@ -103,7 +104,7 @@ bool initial_guess::sad::setup(OrbitalVector &Phi, double prec, double screen, c initial_guess::sad::project_atomic_densities(prec, rho_j, nucs, screen); // Compute XC density - Density &rho_xc = XC.getDensity(DensityType::Total); + Density &rho_xc = XC_.getDensity(DensityType::Total); mrcpp::deep_copy(rho_xc, rho_j); if (plevel == 1) mrcpp::print::time(1, "Projecting GTO density", t_lap); @@ -151,6 +152,7 @@ bool initial_guess::sad::setup(OrbitalVector &Phi, double prec, double screen, c print_utils::text(0, "Screening ", print_utils::dbl_to_str(screen, 5, true) + " StdDev"); print_utils::text(0, "Restricted ", (restricted) ? "True" : "False"); print_utils::text(0, "Functional ", "LDA (SVWN5)"); + print_utils::text(0, "XC Library ", (mrdft::Factory::libxc) ? "LibXC" : "XCFun"); print_utils::text(0, "AO basis ", "3-21G"); mrcpp::print::separator(0, '~', 2); @@ -167,8 +169,8 @@ bool initial_guess::sad::setup(OrbitalVector &Phi, double prec, double screen, c MomentumOperator p(D_p); NuclearOperator V_nuc(nucs, prec); CoulombOperator J(P_p); - XCOperator XC(mrdft_p); - RankZeroOperator V = V_nuc + J + XC; + XCOperator XC_(mrdft_p); + RankZeroOperator V = V_nuc + J + XC_; auto plevel = Printer::getPrintLevel(); if (plevel == 1) mrcpp::print::header(1, "SAD Initial Guess"); @@ -180,7 +182,7 @@ bool initial_guess::sad::setup(OrbitalVector &Phi, double prec, double screen, c initial_guess::sad::project_atomic_densities(prec, rho_j, nucs, screen); // Compute XC density - Density &rho_xc = XC.getDensity(DensityType::Total); + Density &rho_xc = XC_.getDensity(DensityType::Total); mrcpp::deep_copy(rho_xc, rho_j); if (plevel == 1) mrcpp::print::time(1, "Projecting GTO density", t_lap); diff --git a/src/mrdft/Factory.cpp b/src/mrdft/Factory.cpp index fb7ee1189..0e3b35fb1 100644 --- a/src/mrdft/Factory.cpp +++ b/src/mrdft/Factory.cpp @@ -42,23 +42,144 @@ Factory::Factory(const mrcpp::MultiResolutionAnalysis<3> &MRA) : mra(MRA) , xcfun_p(xcfun_new(), xcfun_delete) {} -/** @brief Build a MRDFT object from the currently defined parameters */ +bool Factory::libxc; + +void MapFuncName(std::string name, std::vector &ids, std::vector &coefs) { + // ensure name is upper case + std::transform(name.begin(), name.end(), name.begin(), [](unsigned char c) { return std::toupper(c); }); + + if (name == "PBE0") { + ids = {XC_HYB_GGA_XC_PBEH}; + coefs = {1.0}; + return; + } else if (name == "PBE") { + ids = {XC_GGA_X_PBE, XC_GGA_C_PBE}; + coefs = {1.0, 1.0}; + return; + } else if (name == "SLATERX") { + ids = {XC_LDA_X}; + coefs = {1.0}; + return; + } else if (name == "BECKEX") { + ids = {XC_GGA_X_B88}; + coefs = {1.0}; + return; + } else if (name == "VWN5C") { + ids = {XC_LDA_C_VWN}; + coefs = {1.0}; + return; + } else if (name == "SVWN5") { + ids = {XC_LDA_C_VWN, XC_LDA_X}; + coefs = {1.0, 1.0}; + return; + } else if (name == "B3P86") { + ids = {XC_HYB_GGA_XC_B3P86}; + coefs = {1.0}; + return; + } else if (name == "BPW91") { + ids = {XC_GGA_X_B88, XC_GGA_C_PW91}; + coefs = {1.0, 1.0}; + return; + } else if (name == "B3LYP") { + // Keep as b3lyp5 for now to be consistent with xcfun + // TODO: change the definition of b3lyp in mrchem to not be b3lyp5 + ids = {XC_HYB_GGA_XC_B3LYP5}; + // ids = {XC_HYB_GGA_XC_B3LYP}; + coefs = {1.0}; + return; + } else if (name == "B3LYP5") { + ids = {XC_HYB_GGA_XC_B3LYP5}; + coefs = {1.0}; + return; + } else { + // Check if Libxc has this functional + int number = xc_functional_get_number(name.c_str()); + if (number == -1) { MSG_ABORT(name + " is not a known shorthand in MRChem nor a functional in Libxc!\n"); } + + ids = {number}; + coefs = {1.0}; + return; + } +} + + +void Factory::setFunctional(const std::string &name, double c) { + setLibxc(libxc); // should probably be where setFunctional is called + + if (libxc) { + std::vector ids; + std::vector coefs; + + MapFuncName(name, ids, coefs); + xc_func_type libxc_obj; + for (size_t i = 0; i < ids.size(); i++) { + auto return_code = xc_func_init(&libxc_obj, ids[i], spin ? XC_POLARIZED : XC_UNPOLARIZED); + if (return_code != 0) { + std::ostringstream oss; + oss << "Functional id = " << ids[i] << " not found in the employed version of Libxc.\n"; + MSG_ABORT(oss.str()); + } + xc_func_set_dens_threshold(&libxc_obj, cutoff); + libxc_objects.push_back(libxc_obj); + libxc_coefs.push_back(c * coefs[i]); + } + + } else { + xcfun_set(xcfun_p.get(), name.c_str(), c); + } +} + std::unique_ptr Factory::build() { // Init DFT grid auto grid_p = std::make_unique(mra); + setLibxc(libxc); + + // Init XCFun or Libxc + bool gga = false; + if (libxc) { + for (const auto &f : libxc_objects) { + + switch (f.info->family) { + case XC_FAMILY_LDA: +#ifdef XC_FAMILY_HYB_GGA + case XC_FAMILY_HYB_LDA: +#endif + gga = false; + break; + + case XC_FAMILY_GGA: +#ifdef XC_FAMILY_HYB_GGA + case XC_FAMILY_HYB_GGA: +#endif + gga = true; + break; + + case XC_FAMILY_MGGA: + case XC_FAMILY_HYB_MGGA: + MSG_ABORT("Meta-GGA functionals are not supported in MRChem.!\n"); + + default: + MSG_ABORT("Libxc functional family not handled in MRChem!\n"); + } + + // Check if functional is range separated + if ((f.info->flags & XC_FLAGS_HYB_CAM) || (f.info->flags & XC_FLAGS_HYB_LC)) MSG_ABORT("Coulomb attenuated functionals not supported in MRChem!\n"); + if ((f.info->flags & XC_FLAGS_HYB_CAMY) || (f.info->flags & XC_FLAGS_HYB_LCY)) MSG_ABORT("Yukawa attenuated functionals not supported in MRChem!\n"); + } + } else { + gga = xcfun_is_gga(xcfun_p.get()); + unsigned int mode = 1; //!< only partial derivative mode implemented + unsigned int func_type = (gga) ? 1 : 0; //!< only LDA and GGA supported for now + unsigned int dens_type = 1 + spin; //!< only n (dens_type = 1) or alpha & beta (denst_type = 2) supported now. + unsigned int laplacian = 0; //!< no laplacian + unsigned int kinetic = 0; //!< no kinetic energy density + unsigned int current = 0; //!< no current density + unsigned int exp_derivative = not(gamma); //!< use gamma or explicit derivatives + if (not(gga)) exp_derivative = 0; //!< fall back to gamma-type derivatives if LDA + xcfun_user_eval_setup(xcfun_p.get(), order, func_type, dens_type, mode, laplacian, kinetic, current, exp_derivative); + } - // Init XCFun - bool gga = xcfun_is_gga(xcfun_p.get()); - bool lda = not(gga); - unsigned int mode = 1; //!< only partial derivative mode implemented - unsigned int func_type = (gga) ? 1 : 0; //!< only LDA and GGA supported for now - unsigned int dens_type = 1 + spin; //!< only n (dens_type = 1) or alpha & beta (denst_type = 2) supported now. - unsigned int laplacian = 0; //!< no laplacian - unsigned int kinetic = 0; //!< no kinetic energy density - unsigned int current = 0; //!< no current density - unsigned int exp_derivative = not(gamma); //!< use gamma or explicit derivatives - if (not(gga)) exp_derivative = 0; //!< fall back to gamma-type derivatives if LDA - xcfun_user_eval_setup(xcfun_p.get(), order, func_type, dens_type, mode, laplacian, kinetic, current, exp_derivative); + bool lda = not gga; // Init MW derivative if (gga) { @@ -71,12 +192,15 @@ std::unique_ptr Factory::build() { std::unique_ptr func_p{nullptr}; if (spin) { if (gga) func_p = std::make_unique(order, xcfun_p, diff_p); - if (lda) func_p = std::make_unique(order, xcfun_p); + else if (lda) func_p = std::make_unique(order, xcfun_p); + else MSG_ABORT("Case not handled"); } else { if (gga) func_p = std::make_unique(order, xcfun_p, diff_p); - if (lda) func_p = std::make_unique(order, xcfun_p); + else if (lda) func_p = std::make_unique(order, xcfun_p); + else MSG_ABORT("Case not handled"); } if (func_p == nullptr) MSG_ABORT("Invalid functional type"); + if (libxc) { func_p->set_libxc_functional_object(libxc_objects, libxc_coefs); } diff_p = std::make_unique>(mra, 0.0, 0.0); func_p->setDerivOp(diff_p); func_p->setLogGradient(log_grad); diff --git a/src/mrdft/Factory.h b/src/mrdft/Factory.h index ce9caf556..a2bec6a9c 100644 --- a/src/mrdft/Factory.h +++ b/src/mrdft/Factory.h @@ -27,37 +27,108 @@ #include #include +#include +#include #include "MRDFT.h" namespace mrdft { +/** + * @class Factory + * @brief Building class for MRDFT objects + * @details Manages the different xc libraries, mapping of functional names and builds the functional objects + * with the required parameters to initialize a DFT calculation + */ class Factory final { public: + /** + * @brief Construct a new Factory object. Initializes the MRA reference and creates the XCFun handle + * @param[in] MRA The Multi-Resolution Analysis object providing the grid + * and basis functions for the calculation + */ Factory(const mrcpp::MultiResolutionAnalysis<3> &MRA); - ~Factory() = default; - void setSpin(bool s) { spin = s; } - void setOrder(int k) { order = k; } - void setUseGamma(bool g) { gamma = g; } - void setLogGradient(bool lg) { log_grad = lg; } - void setDensityCutoff(double c) { cutoff = c; } - void setDerivative(const std::string &n) { diff_s = n; } - void setFunctional(const std::string &n, double c = 1.0) { xcfun_set(xcfun_p.get(), n.c_str(), c); } + ~Factory() = default; ///< @brief Default destructor + + /* + * Setters + */ + void setSpin(bool s) { spin = s; } ///< Set spin polarization (true for unrestricted/spin-polarized) */ + void setOrder(int k) { order = k; } ///< Set the polynomial order for the MRA basis + void setUseGamma(bool g) { gamma = g; } ///< Toggle between gamma-type and explicit derivatives + void setLogGradient(bool lg) { log_grad = lg; } ///< Toggle the use of logarithmic gradients + void setDensityCutoff(double c) { cutoff = c; } ///< Set the threshold for neglecting low-density regions + void setLibxc(bool libxc_) {libxc = libxc_; } ///< Toggle between Libxc (true) and XCFun (false) backends + void setDerivative(const std::string &n) { diff_s = n; } ///< Set derivative operator type (e.g., "bspline", "abgv_00") + + /** + * @brief Configures the xc functional + * + * @param[in] name The name of the xc functional (e.g., "PBE", "B3LYP") + * @param[in] c A global scaling coefficient applied to the functional. + * @throws MSG_ABORT If a mapped Libxc ID is incompatible with the linked Libxc version + * @note Depending on the chosen library, this method either initializes + * one or more Libxc functional objects using the MapFuncName function or sets + * the functional parameters in the XCFun backend + */ + void setFunctional(const std::string &n, double c = 1.0); + /** + * @brief Build a MRDFT object from the currently defined parameters + * @details Performs the following steps: + * 1. **Grid Initialization**: Creates a multi-resolution grid based on the MRA + * 2. **Library dependent initiation**: If using Libxc, iterates through functional objects + * to ensure they belong to supported families (LDA/GGA, not meta-GGA or range separated). + * If using XCFun, sets evaluation parameters, mode and order + * 3. **Operator Selection**: Assigns numerical derivative operators (BSpline or ABGV) + * required for GGAs + * 4. **Functional Instantiation**: Selects the appropriate concrete implementation + * (SpinLDA, SpinGGA, LDA, or GGA) based on spin and gradient requirements. + * 5. **State Sync**: Passes functional objects, density cutoffs, and derivative + * schemes to the functional + * @return std::unique_ptr A pointer to the assembled Multi-Resolution DFT object. + * @throws MSG_ABORT If unsupported functional families are detected in the Libxc case + * (eg. meta-GGAs and range separated functionals) + */ std::unique_ptr build(); + static bool libxc; ///< @brief Flag indicating if Libxc is active (True if "DFT {xc_library = libxc}" in input file) + +private: private: - int order{1}; - bool spin{false}; - bool gamma{false}; - bool log_grad{false}; - double cutoff{-1.0}; - std::string diff_s{"abgv_00"}; - const mrcpp::MultiResolutionAnalysis<3> mra; - - XC_p xcfun_p; - std::unique_ptr> diff_p; + int order{1}; ///< Polynomial order of the Multi-Resolution Analysis (MRA) basis + bool spin{false}; ///< If true, perform unrestricted calculations + bool gamma{false}; ///< If true, use gamma-type derivatives (gradient squared) instead of explicit components + bool log_grad{false}; ///< Toggle for using logarithmic gradient transformations + double cutoff{-1.0}; ///< Density threshold; values below this are sat to 0 + std::string diff_s{"abgv_00"}; ///< String identifier for the derivative operator type (e.g., "bspline", "abgv_55") + + const mrcpp::MultiResolutionAnalysis<3> mra; ///< @brief Reference to the 3D Multi-Resolution Analysis grid structure + std::unique_ptr> diff_p; ///< @brief Pointer to the numerical derivative operator used for GGA gradients + XC_p xcfun_p; ///< @brief Pointer to the XCFun library handle + + + std::vector libxc_objects; ///< @brief Vector of initialized Libxc functionals + std::vector libxc_coefs; ///< @brief Vector scaling coefficients for each functional in libxc_objects + + /** + * @brief Maps a functional name string (e.g., "PBE0", "LDA" or "XC_LDA_X", XC_GGA_X_B88) + * to its corresponding Libxc IDs and scaling coefficients + * @note The input `name` is transformed to uppercase internally, making the + * search case-insensitive + * @param[in] name Name of the functional + * @param[in] ids Vector to be populated with the IDs used by Libxc + * @param[in] coefs Vector to be populated with the corresponding scaling coefficients + * @throws MSG_ABORT If the name is not a recognized internal shorthand and + * is not found within the Libxc library + * @example + * std::vector ids; + * std::vector coefs; + * MapFuncName("LDA", ids, coefs); + * // ids: {XC_LDA_C_VWN, XC_LDA_X}, coefs: {1.0, 1.0} + */ + std::vector mapFunctionalName(const std::string &name) const; }; } // namespace mrdft diff --git a/src/mrdft/Functional.cpp b/src/mrdft/Functional.cpp index 35876d6fc..63e51285e 100644 --- a/src/mrdft/Functional.cpp +++ b/src/mrdft/Functional.cpp @@ -24,83 +24,265 @@ */ #include +#include "utils/print_utils.h" +#include +#include "Factory.h" #include "Functional.h" namespace mrdft { -/** @brief Run a collection of grid points through XCFun - * - * Each row corresponds to one grid point. - * - * param[in] inp_data Matrix of input values - * param[out] out_data Matrix of output values - */ -Eigen::MatrixXd Functional::evaluate(Eigen::MatrixXd &inp) const { - int nInp = xcfun_input_length(xcfun.get()); // Input parameters to XCFun - int nOut = xcfun_output_length(xcfun.get()); // Input parameters to XCFun - int nPts = inp.cols(); - if (nInp != inp.rows()) MSG_ABORT("Invalid input"); +void Functional::print_functional_references() const { + // Only run on main thread + if (mrcpp::mpi::world_rank != 0) { + return; + } - Eigen::MatrixXd out = Eigen::MatrixXd::Zero(nOut, nPts); - for (int i = 0; i < nPts; i++) { - bool calc = true; - if (isSpin()) { - if (inp(0, i) < cutoff and inp(1, i) < cutoff) calc = false; - } else { - if (inp(0, i) < cutoff) calc = false; + auto pwidth = mrcpp::Printer::getWidth(); + auto txt_width = 50; + auto pre_spaces = (pwidth - 6 - txt_width) / 2; + auto post_spaces = pwidth - 6 - txt_width - pre_spaces; + std::string pre_str = std::string(3, '*') + std::string(pre_spaces, ' '); + std::string post_str = std::string(post_spaces, ' ') + std::string(3, '*'); + + mrcpp::print::separator(0, ' '); + mrcpp::print::separator(0, ' '); + mrcpp::print::separator(0, '*'); + println(0, pre_str << " " << post_str); + println(0, pre_str << " XC Functional " << post_str); + println(0, pre_str << " " << post_str); + mrcpp::print::separator(0, '*', 1); + + // XCFun is used + if (not Factory::libxc) { + printout(0, xcfun_splash()); + mrcpp::print::separator(0, ' '); + mrcpp::print::separator(0, '-', 1); + return; + } + + // Conditional reference printing + auto print_wrap = [&](std::string str, std::size_t txt_width) { + size_t offset = 0; + while (offset + txt_width < str.size()) { + // Is the line already sufficiently short? + size_t newline_pos = str.find('\n', offset); + if (newline_pos < offset + txt_width) { + offset = newline_pos + 1; + continue; + } + + size_t space_pos = str.rfind(' ', offset + txt_width); + // If the string doesn't have a space, or it is too far out, hard insert a newline + if (space_pos == std::string::npos || space_pos - offset > txt_width) { + space_pos = offset + txt_width; + str.insert(space_pos, "\n"); + } else { + str[space_pos] = '\n'; + } + offset = space_pos + 1; + } + std::cout << str; + }; + + std::string str = "Using Libxc (version " + std::string(xc_version_string()) + ") to evaluate density functionals. Libxc is free software. It is " + + "distributed under the Mozilla Public License, version 2.0. For " + "more information, please check the Libxc manual. You should cite\n\n" + xc_reference() + + " DOI: " + xc_reference_doi() + "\n\nwhen " + "reporting the results of your calculation in a scientific article.\n"; + auto libxc_txt_width = 75; + print_wrap(str, libxc_txt_width); + + // Avoid printing the same LibXC functional multiple times + std::set printed_ids; + std::cout << "\nLibxc functionals used in this calculation:\n"; + for (const auto &func : libxc_objects) { + if (func.info == nullptr) continue; // safety + + int id = func.info->number; + if (!printed_ids.insert(id).second) { + // Already printed this ID + continue; + } + + char *name = xc_functional_get_name(id); + std::cout << " - " << name << " (ID " << id << "): " << func.info->name << std::endl; + free(name); + + for (int number = 0; number < XC_MAX_REFERENCES; number++) { + auto reference = xc_func_info_get_references(func.info, number); + if (reference == nullptr) break; + std::cout << " * " << xc_func_reference_get_ref(reference) << "\n DOI:" << xc_func_reference_get_doi(reference) << std::endl; } - // NB: the data is stored colomn major, i.e. two consecutive points of for example energy density, are not consecutive in memory - // That means that we cannot extract the energy density data with out.row(0).data() for example. - if (calc) xcfun_eval(xcfun.get(), inp.col(i).data(), out.col(i).data()); } - return out; } +void Functional::set_libxc_functional_object(std::vector libxc_objects_, std::vector libxc_coefs_) { + libxc_objects = std::move(libxc_objects_); + libxc_coefs = std::move(libxc_coefs_); +} -/** @brief Run a collection of grid points through XCFun - * - * Each column corresponds to one grid point. - * From a performance point of view, (in pre and postprocessing) it is much more - * efficient to have the two consecutive points in two consecutive adresses in memory - * - * param[in] inp_data Matrix of input values - * param[out] out_data Matrix of output values - */ -Eigen::MatrixXd Functional::evaluate_transposed(Eigen::MatrixXd &inp) const { - int nInp = xcfun_input_length(xcfun.get()); // Input parameters to XCFun - int nOut = xcfun_output_length(xcfun.get()); // Input parameters to XCFun - int nPts = inp.rows(); - if (nInp != inp.cols()) MSG_ABORT("Invalid input"); - - Eigen::MatrixXd out = Eigen::MatrixXd::Zero(nPts, nOut); - Eigen::VectorXd inp_row = Eigen::VectorXd::Zero(nInp); - Eigen::VectorXd out_row = Eigen::VectorXd::Zero(nOut); - for (int i = 0; i < nPts; i++) { - bool calc = true; - if (isSpin()) { - if (inp(i, 0) < cutoff and inp(i, 1) < cutoff) calc = false; - } else { - if (inp(i, 0) < cutoff) calc = false; +double Functional::amountEXX() const { + double exx = 0.0; + if (Factory::libxc) { + for (std::size_t i = 0; i < libxc_objects.size(); ++i) { + const xc_func_type &f = libxc_objects[i]; + double frac = xc_hyb_exx_coef(&f); + exx += libxc_coefs[i] * frac; } - if (calc) { - inp_row = inp.row(i); - xcfun_eval(xcfun.get(), inp_row.data(), out_row.data()); - out.row(i) = out_row; + } else { + xcfun_get(xcfun.get(), "exx", &exx); + } + return exx; +} + +void Functional::evaluate_data(const Eigen::MatrixXd &inp, Eigen::MatrixXd &out) const { + int nInp = numIn(); + int nOut = numOut(); + int nPts = inp.cols(); + if (nInp != inp.rows()) { + std::ostringstream oss; + oss << "Invalid input: expected matrix with " << nInp << " rows, got " << inp.rows() << "!\n"; + MSG_ABORT(oss.str()); + } + if (nOut != out.rows()) { + std::ostringstream oss; + oss << "Invalid output: expected matrix with " << nOut << " rows, got " << out.rows() << "!\n"; + MSG_ABORT(oss.str()); + } + out.setZero(); + + if (Factory::libxc) { + Eigen::MatrixXd exc, vxc, sxc, sigma; + for (size_t i = 0; i < libxc_objects.size(); i++) { + switch (libxc_objects[i].info->family) { + case XC_FAMILY_LDA: + case XC_FAMILY_HYB_LDA: + if (isSpin()) { + Eigen::MatrixXd rho = Eigen::MatrixXd::Zero(2, nPts); + exc = Eigen::MatrixXd::Zero(1, nPts); + vxc = Eigen::MatrixXd::Zero(2, nPts); + for (size_t j = 0; j < nPts; j++) { + // alpha_1, beta_1, alpha_2, beta_2, .. + rho(0, j) = inp(0, j); + rho(1, j) = inp(1, j); + } + xc_lda_exc_vxc(&libxc_objects[i], nPts, rho.data(), exc.data(), vxc.data()); + for (size_t j = 0; j < nPts; ++j) { + // xcfun calculates actual energy density while libxc calculates + // energy density per electron density + + // rho = rho_alpha + rho_beta + out(0, j) += exc(0, j) * libxc_coefs[i] * (inp(0, j) + inp(1, j)); + out(1, j) += vxc(0, j) * libxc_coefs[i]; + out(2, j) += vxc(1, j) * libxc_coefs[i]; + } + } else { + exc = Eigen::MatrixXd::Zero(1, nPts); + vxc = Eigen::MatrixXd::Zero(1, nPts); + xc_lda_exc_vxc(&libxc_objects[i], nPts, inp.data(), exc.data(), vxc.data()); + for (size_t j = 0; j < nPts; ++j) { + // xcfun calculates actual energy density while libxc calculates + // energy density per electron density + out(0, j) += exc(0, j) * libxc_coefs[i] * inp(0, j); + out(1, j) += vxc(0, j) * libxc_coefs[i]; + } + } + break; + case XC_FAMILY_GGA: + case XC_FAMILY_HYB_GGA: + if (isSpin()) { + Eigen::MatrixXd rho = Eigen::MatrixXd::Zero(2, nPts); + exc = Eigen::MatrixXd::Zero(1, nPts); + vxc = Eigen::MatrixXd::Zero(2, nPts); + sxc = Eigen::MatrixXd::Zero(3, nPts); + sigma = Eigen::MatrixXd::Zero(3, nPts); + for (size_t j = 0; j < nPts; j++) { + // alpha_1, beta_1, alpha_2, beta_2, .. + rho(0, j) = inp(0, j); + rho(1, j) = inp(1, j); + } + for (size_t j = 0; j < nPts; j++) { + // clang-format off + // Libxc takes in reduced gradients: up-up, up-down, down-down + sigma(0, j) = inp(2, j) * inp(2, j) + inp(3, j) * inp(3, j) + inp(4, j) * inp(4, j); + sigma(1, j) = inp(2, j) * inp(5, j) + inp(3, j) * inp(6, j) + inp(4, j) * inp(7, j); + sigma(2, j) = inp(5, j) * inp(5, j) + inp(6, j) * inp(6, j) + inp(7, j) * inp(7, j); + } + xc_gga_exc_vxc(&libxc_objects[i], nPts, rho.data(), sigma.data(), exc.data(), vxc.data(), sxc.data()); + + for (size_t j = 0; j < nPts; ++j) { + // clang-format off + // xcfun calculates energy density per volume while libxc calculates + // energy density per electron, so we multiply by the density here + out(0, j) += exc(0, j) * libxc_coefs[i] * (inp(0, j) + inp(1, j)); + out(1, j) += vxc(0, j) * libxc_coefs[i]; + out(2, j) += vxc(1, j) * libxc_coefs[i]; + + // alpha_i, coef * ( 2 * vaa * grad_a_i + vab * grad_b_i ), i = x, y, z + out(3, j) += libxc_coefs[i] * ( 2 * sxc(0, j) * inp(2, j) + sxc(1, j) * inp(5, j) ); + out(4, j) += libxc_coefs[i] * ( 2 * sxc(0, j) * inp(3, j) + sxc(1, j) * inp(6, j) ); + out(5, j) += libxc_coefs[i] * ( 2 * sxc(0, j) * inp(4, j) + sxc(1, j) * inp(7, j) ); + // beta_i, coef * ( 2 * vbb * grad_b_i + vab * grad_a_i ), i = x, y, z + out(6, j) += libxc_coefs[i] * ( 2 * sxc(2, j) * inp(5, j) + sxc(1, j) * inp(2, j) ); + out(7, j) += libxc_coefs[i] * ( 2 * sxc(2, j) * inp(6, j) + sxc(1, j) * inp(3, j) ); + out(8, j) += libxc_coefs[i] * ( 2 * sxc(2, j) * inp(7, j) + sxc(1, j) * inp(4, j) ); + // clang-format on + } + } else { + Eigen::MatrixXd rho = inp.row(0).transpose(); + exc = Eigen::MatrixXd::Zero(1, nPts); + vxc = Eigen::MatrixXd::Zero(1, nPts); + sxc = Eigen::MatrixXd::Zero(1, nPts); + sigma = Eigen::MatrixXd::Zero(1, nPts); + for (size_t j = 0; j < nPts; j++) { sigma(j) = inp(1, j) * inp(1, j) + inp(2, j) * inp(2, j) + inp(3, j) * inp(3, j); } + xc_gga_exc_vxc(&libxc_objects[i], nPts, rho.data(), sigma.data(), exc.data(), vxc.data(), sxc.data()); + + for (size_t j = 0; j < nPts; ++j) { + // xcfun calculates energy density per volume while libxc calculates + // energy density per electron, so we multiply by the density here + out(0, j) += exc(0, j) * libxc_coefs[i] * inp(0, j); + out(1, j) += vxc(0, j) * libxc_coefs[i]; + out(2, j) += 2 * sxc(0, j) * inp(1, j) * libxc_coefs[i]; + out(3, j) += 2 * sxc(0, j) * inp(2, j) * libxc_coefs[i]; + out(4, j) += 2 * sxc(0, j) * inp(3, j) * libxc_coefs[i]; + } + } + break; + default: + break; + } + } + } else { + if (nInp != xcfun_input_length(xcfun.get()) or nOut != xcfun_output_length(xcfun.get())) { throw std::logic_error("Dimension mismatch!\n"); } + + for (int i = 0; i < nPts; i++) { + if (isSpin()) { + if (inp(0, i) < cutoff and inp(1, i) < cutoff) continue; + } else { + if (inp(0, i) < cutoff) continue; + } + xcfun_eval(xcfun.get(), inp.col(i).data(), out.col(i).data()); } } +} + +Eigen::MatrixXd Functional::evaluate(Eigen::MatrixXd &inp) const { + int nOut = numOut(); + int nPts = inp.cols(); + + Eigen::MatrixXd out = Eigen::MatrixXd::Zero(nOut, nPts); + evaluate_data(inp, out); return out; } +Eigen::MatrixXd Functional::evaluate_transposed(Eigen::MatrixXd &inp) const { + // NB: the data is stored colomn major, i.e. two consecutive points of for example energy density, are not consecutive in memory + // That means that we cannot extract the energy density data with out.row(0).data() for example. + Eigen::MatrixXd inp_trans(inp.transpose()); + Eigen::MatrixXd out_trans(numOut(), inp.rows()); + evaluate_data(inp_trans, out_trans); + return out_trans.transpose(); +} -/** @brief Contract a collection of grid points - * - * Each row corresponds to one grid point. - * - * param[in] xc_data Matrix of functional partial derivative values - * param[in] d_data Matrix of density input values - * param[out] out_data Matrix of contracted output values - */ Eigen::MatrixXd Functional::contract(Eigen::MatrixXd &xc_data, Eigen::MatrixXd &d_data) const { auto nPts = xc_data.cols(); auto nFcs = getCtrOutputLength(); @@ -125,14 +307,6 @@ Eigen::MatrixXd Functional::contract(Eigen::MatrixXd &xc_data, Eigen::MatrixXd & return out_data; } -/** @brief Contract a collection of grid points - * - * Each column corresponds to one set of grid points. - * - * param[in] xc_data Matrix of functional partial derivative values - * param[in] d_data Matrix of density input values - * param[out] out_data Matrix of contracted output values - */ Eigen::MatrixXd Functional::contract_transposed(Eigen::MatrixXd &xc_data, Eigen::MatrixXd &d_data) const { auto nPts = xc_data.rows(); auto nFcs = getCtrOutputLength(); @@ -156,56 +330,6 @@ Eigen::MatrixXd Functional::contract_transposed(Eigen::MatrixXd &xc_data, Eigen: return out_data; } - -/** @brief Evaluates XC functional and derivatives for a given NodeIndex - * - * The electronic densities (total/alpha/beta) are given as input. - * The values of the zero order densities and their gradient are sent to xcfun. - * The output of xcfun must then be combined ("contract") with the gradients - * of the higher order densities. - * - * XCFunctional output (with k=1 and explicit derivatives): - * - * LDA: \f$ \left(F_{xc}, \frac{\partial F_{xc}}{\partial \rho}\right) \f$ - * - * GGA: \f$ \left(F_{xc}, - * \frac{\partial F_{xc}}{\partial \rho}, - * \frac{\partial F_{xc}}{\partial \rho_x}, - * \frac{\partial F_{xc}}{\partial \rho_y}, - * \frac{\partial F_{xc}}{\partial \rho_z}\right) \f$ - * - * Spin LDA: \f$ \left(F_{xc}, \frac{\partial F_{xc}}{\partial \rho^\alpha}, - * \frac{\partial F_{xc}}{\partial \rho^\beta}\right) \f$ - * - * Spin GGA: \f$ \left(F_{xc}, - * \frac{\partial F_{xc}}{\partial \rho^\alpha}, - * \frac{\partial F_{xc}}{\partial \rho^\beta}, - * \frac{\partial F_{xc}}{\partial \rho_x^\alpha}, - * \frac{\partial F_{xc}}{\partial \rho_y^\alpha}, - * \frac{\partial F_{xc}}{\partial \rho_z^\alpha}, - * \frac{\partial F_{xc}}{\partial \rho_x^\beta}, - * \frac{\partial F_{xc}}{\partial \rho_y^\beta}, - * \frac{\partial F_{xc}}{\partial \rho_z^\beta} - * \right) \f$ - * - * XCFunctional output (with k=1 and gamma-type derivatives): - * - * GGA: \f$ \left(F_{xc}, - * \frac{\partial F_{xc}}{\partial \rho}, - * \frac{\partial F_{xc}}{\partial \gamma} \f$ - * - * Spin GGA: \f$ \left(F_{xc}, - * \frac{\partial F_{xc}}{\partial \rho^\alpha}, - * \frac{\partial F_{xc}}{\partial \rho^\beta }, - * \frac{\partial F_{xc}}{\partial \gamma^{\alpha \alpha}}, - * \frac{\partial F_{xc}}{\partial \gamma^{\alpha \beta }}, - * \frac{\partial F_{xc}}{\partial \gamma^{\beta \beta }} - * \right) \f$ - * - * param[in] inp Input values - * param[out] xcNodes Output values - * - */ void Functional::makepot(mrcpp::FunctionTreeVector<3> &inp, std::vector *> xcNodes) const { if (this->log_grad){ MSG_ERROR("log_grad not implemented"); @@ -248,7 +372,7 @@ void Functional::makepot(mrcpp::FunctionTreeVector<3> &inp, std::vector - -#include -#include #include #include #include +#include +#include namespace mrdft { using XC_p = std::unique_ptr; +/** + * @class Functional + * @brief Abstract base class for Exchange-Correlation functionals + * @details Provides the interface for evaluating XC potentials + * on the Multi-Resolution Analysis (MRA) grid using either Libxc or XCFun + */ class Functional { public: + /** + * @brief Constructor for the Functional base class + * @param[in] k The polynomial order for the MRA basis + * @param[in] f The XCFun handle (ownership is transferred) + */ Functional(int k, XC_p &f) : order(k) , xcfun(std::move(f)) {} virtual ~Functional() = default; + /** @brief Evaluates XC functional and derivatives for a given NodeIndex + * + * The electronic densities (total/alpha/beta) are given as input. + * The values of the zero order densities and their gradient are sent to xcfun. + * The output of xcfun must then be combined ("contract") with the gradients + * of the higher order densities. + * + * XCFunctional output (with k=1 and explicit derivatives): + * + * LDA: \f$ \left(F_{xc}, \frac{\partial F_{xc}}{\partial \rho}\right) \f$ + * + * GGA: \f$ \left(F_{xc}, + * \frac{\partial F_{xc}}{\partial \rho}, + * \frac{\partial F_{xc}}{\partial \rho_x}, + * \frac{\partial F_{xc}}{\partial \rho_y}, + * \frac{\partial F_{xc}}{\partial \rho_z}\right) \f$ + * + * Spin LDA: \f$ \left(F_{xc}, \frac{\partial F_{xc}}{\partial \rho^\alpha}, + * \frac{\partial F_{xc}}{\partial \rho^\beta}\right) \f$ + * + * Spin GGA: \f$ \left(F_{xc}, + * \frac{\partial F_{xc}}{\partial \rho^\alpha}, + * \frac{\partial F_{xc}}{\partial \rho^\beta}, + * \frac{\partial F_{xc}}{\partial \rho_x^\alpha}, + * \frac{\partial F_{xc}}{\partial \rho_y^\alpha}, + * \frac{\partial F_{xc}}{\partial \rho_z^\alpha}, + * \frac{\partial F_{xc}}{\partial \rho_x^\beta}, + * \frac{\partial F_{xc}}{\partial \rho_y^\beta}, + * \frac{\partial F_{xc}}{\partial \rho_z^\beta} + * \right) \f$ + * + * XCFunctional output (with k=1 and gamma-type derivatives): + * + * GGA: \f$ \left(F_{xc}, + * \frac{\partial F_{xc}}{\partial \rho}, + * \frac{\partial F_{xc}}{\partial \gamma} \f$ + * + * Spin GGA: \f$ \left(F_{xc}, + * \frac{\partial F_{xc}}{\partial \rho^\alpha}, + * \frac{\partial F_{xc}}{\partial \rho^\beta }, + * \frac{\partial F_{xc}}{\partial \gamma^{\alpha \alpha}}, + * \frac{\partial F_{xc}}{\partial \gamma^{\alpha \beta }}, + * \frac{\partial F_{xc}}{\partial \gamma^{\beta \beta }} + * \right) \f$ + * + * @param[in] inp Input values + * @param[out] xcNodes Output values + * + */ void makepot(mrcpp::FunctionTreeVector<3> &inp, std::vector *> xcNodes) const; - void setLogGradient(bool log) { log_grad = log; } - void setDensityCutoff(double cut) { cutoff = cut; } - void setDerivOp(std::unique_ptr> &d) { derivOp = std::move(d); } - - virtual bool isSpin() const = 0; - bool isLDA() const { return (not(isGGA() or isMetaGGA())); } - bool isGGA() const { return xcfun_is_gga(xcfun.get()); } - bool isMetaGGA() const { return xcfun_is_metagga(xcfun.get()); } - bool isHybrid() const { return (std::abs(amountEXX()) > 1.0e-10); } - double amountEXX() const { - double exx = 0.0; - xcfun_get(xcfun.get(), "exx", &exx); - return exx; - } - double XCenergy = 0.0; - + /** + * Setters + */ + void setLogGradient(bool log) { log_grad = log; } ///< @brief Set whether to use logarithmic gradient transformations + void setDensityCutoff(double cut) { cutoff = cut; } ///< @brief Set the density threshold below which density is set to 0 + void setDerivOp(std::unique_ptr> &d) { derivOp = std::move(d); } ///< @brief Set the numerical derivative operator for gradient-based functionals + + /** + * Functional type querying + */ + bool isLDA() const { return not (isGGA() or isMetaGGA()); } ///< @return True if functional is LDA type (not a GGA or meta-GGA) + bool isHybrid() const { return (std::abs(amountEXX()) > 1.0e-10); } ///< @return True if functional is a hybrid (includes exact exchange) + virtual bool isSpin() const = 0; ///< @brief Returns True if the functional object is spin-polarized + virtual bool isGGA() const = 0; ///< @brief Returns True if the functional is a GGA + virtual bool isMetaGGA() const = 0; ///< @brief Returns True if the functional is a Meta-GGA + + virtual int numIn() const = 0; ///< Fetches number of variables in the input matrix + virtual int numOut() const = 0; ///< Fetches number of variables in the output matrix + + /** + * @brief Fetches the amount of exact exchange needed for a given functional + * @return The total fraction of exx to be added to the functional + */ + double amountEXX() const; + double XCenergy = 0.0; ///< @brief Stores calculated xc energy for the current state + + /** + * @brief Evaluates the functional on a set of grid points + * @param[in] inp Matrix of input values (density, gradient,...) + * @param[out] out out_data Matrix of output values (energy, potential, ...) + * @details Each column corresponds to one grid point + * From a performance point of view, (in pre and postprocessing) it is much more + * efficient to have the two consecutive points in two consecutive adresses in memory + */ Eigen::MatrixXd evaluate(Eigen::MatrixXd &inp) const; + + /** + * @brief Evaluates the functional on a set of grid points. Transposed version of Functional::evaluate() + * @param[in] inp Matrix of input values (density, gradient,...) + * @param[out] out_trans Matrix of output values (energy, potential, ...) + * @details Each row corresponds to one grid point + */ Eigen::MatrixXd evaluate_transposed(Eigen::MatrixXd &inp) const; + + bool libxc; ///< @brief Flag indicating if Libxc is active (True if "DFT {xc_library = libxc}" in input file) + std::vector libxc_objects; ///< @brief Vector of initialized Libxc functionals + std::vector libxc_coefs; ///< @brief Vector scaling coefficients for each functional in libxc_objects + + /** + * @brief Prints the splash screens, version info, and references for the + * active xc libraries and functionals + * @details If Libxc is used, it iterates through + * all initialized functional objects to print their specific DOIs + */ + void print_functional_references() const; + + /** + * @brief Transfers ownership of Libxc functional objects and their scaling + * coefficients to the Functional instance + * @param[in] libxc_objects_ Vector of initialized Libxc functionals + * @param[in] libxc_coefs_ Vector of corresponding weights of the initialized Libxc functionals + */ + void set_libxc_functional_object(std::vector libxc_objects_, std::vector libxc_coefs_); + friend class MRDFT; protected: - const int order; - bool log_grad{false}; - double cutoff{-1.0}; - Eigen::VectorXi d_mask; - Eigen::MatrixXi xc_mask; - XC_p xcfun; - std::unique_ptr> derivOp{nullptr}; - - int getXCInputLength() const { return xcfun_input_length(xcfun.get()); } - int getXCOutputLength() const { return xcfun_output_length(xcfun.get()); } - virtual int getCtrInputLength() const = 0; - virtual int getCtrOutputLength() const = 0; - + const int order; ///< @brief Order of functional derivatives. Eg. 0 (energy), 1 (potential), 2 (pot. gradient), etc. + bool log_grad{false}; ///< @brief Toggle for logarithmic gradient + double cutoff{-1.0}; ///< @brief Density threshold + Eigen::VectorXi d_mask; ///< @brief density and derivative(s) mask vector for response calculations + Eigen::MatrixXi xc_mask; ///< @brief functional and derivative(s) mask vector for response calculations + XC_p xcfun; ///< @brief XCFun library handle + std::unique_ptr> derivOp{nullptr}; ///< @brief Operator used to compute gradients + + /** + * @brief Run a collection of grid points through Libxc or XCFun + * @param[in] inp Matrix of input values, where each row is one grid point + * @param[out] out Matrix of output values + */ + void evaluate_data(const Eigen::MatrixXd & inp, Eigen::MatrixXd &out) const; + + /** + * @brief Contracts a collection of grid points + * @details This is used to implement the chain rule for functionals involving + * gradients or when calculating higher-order properties + * @param[in] xc_data xc_data Matrix of functional partial derivative values + * @param[in] d_data d_data Matrix of density input values + * @param[out] out_data Matrix of contracted output values + * @details Each row corresponds to one grid point + */ Eigen::MatrixXd contract(Eigen::MatrixXd &xc_data, Eigen::MatrixXd &d_data) const; + /** + * @brief Contracts a collection of grid points. Transposed version of Functional::contract + * @details This is used to implement the chain rule for functionals involving + * gradients or when calculating higher-order properties + * @param[in] xc_data xc_data Matrix of functional partial derivative values + * @param[in] d_data d_data Matrix of density input values + * @param[out] out_data Matrix of contracted output values + * @details Each column corresponds to one grid point + */ Eigen::MatrixXd contract_transposed(Eigen::MatrixXd &xc_data, Eigen::MatrixXd &d_data) const; - virtual void clear() = 0; - virtual mrcpp::FunctionTreeVector<3> setupXCInput() = 0; - virtual mrcpp::FunctionTreeVector<3> setupCtrInput() = 0; + virtual int getCtrInputLength() const = 0; ///< @brief Expected number of input components for the contraction step + virtual int getCtrOutputLength() const = 0; ///< @brief Expected number of output components for the contraction step + virtual void clear() = 0; ///< @brief Clears internal functions + virtual mrcpp::FunctionTreeVector<3> setupXCInput() = 0; ///< @brief Configures input for evaluation + virtual mrcpp::FunctionTreeVector<3> setupCtrInput() = 0; ///< @brief Configures input for contraction + virtual void preprocess(mrcpp::FunctionTreeVector<3> &inp) = 0; ///< @brief Collects input functions for evaluation + virtual mrcpp::FunctionTreeVector<3> postprocess(mrcpp::FunctionTreeVector<3> &inp) = 0; ///< @brief Computes final output functions - virtual void preprocess(mrcpp::FunctionTreeVector<3> &inp) = 0; - virtual mrcpp::FunctionTreeVector<3> postprocess(mrcpp::FunctionTreeVector<3> &inp) = 0; }; } // namespace mrdft diff --git a/src/mrdft/GGA.h b/src/mrdft/GGA.h index 3ba7ec89c..4379b75d7 100644 --- a/src/mrdft/GGA.h +++ b/src/mrdft/GGA.h @@ -28,6 +28,7 @@ #include #include "Functional.h" +#include "Factory.h" // only to call Factory::libxc namespace mrdft { @@ -37,6 +38,10 @@ class GGA final : public Functional { ~GGA() override = default; bool isSpin() const override { return false; } + bool isGGA() const override { return true; } + bool isMetaGGA() const override { return false; } + int numIn() const override { return 4; } + int numOut() const override { if (Factory::libxc) {return 5;} else {return xcfun_output_length(xcfun.get());} } private: std::unique_ptr> derivative{nullptr}; diff --git a/src/mrdft/LDA.h b/src/mrdft/LDA.h index 88caeb823..f4455f30c 100644 --- a/src/mrdft/LDA.h +++ b/src/mrdft/LDA.h @@ -28,6 +28,7 @@ #include #include "Functional.h" +#include "Factory.h" // only to call Factory::libxc namespace mrdft { @@ -37,6 +38,10 @@ class LDA final : public Functional { ~LDA() override = default; bool isSpin() const override { return false; } + bool isGGA() const override { return false; } + bool isMetaGGA() const override { return false; } + int numIn() const override { return 1; } + int numOut() const override { if (Factory::libxc) {return 2;} else {return xcfun_output_length(xcfun.get());} } private: mrcpp::FunctionTreeVector<3> rho; diff --git a/src/mrdft/SpinGGA.h b/src/mrdft/SpinGGA.h index 9f4613242..34400e9ed 100644 --- a/src/mrdft/SpinGGA.h +++ b/src/mrdft/SpinGGA.h @@ -28,6 +28,7 @@ #include #include "Functional.h" +#include "Factory.h" // only to call Factory::libxc namespace mrdft { @@ -37,6 +38,10 @@ class SpinGGA final : public Functional { ~SpinGGA() override = default; bool isSpin() const override { return true; } + bool isGGA() const override { return true; } + bool isMetaGGA() const override { return false; } + int numIn() const override { return 8; } + int numOut() const override { if (Factory::libxc) {return 9;} else {return xcfun_output_length(xcfun.get());} } private: std::unique_ptr> derivative{nullptr}; diff --git a/src/mrdft/SpinLDA.h b/src/mrdft/SpinLDA.h index 34fca1722..8f17e3971 100644 --- a/src/mrdft/SpinLDA.h +++ b/src/mrdft/SpinLDA.h @@ -28,6 +28,7 @@ #include #include "Functional.h" +#include "Factory.h" // only to call Factory::libxc namespace mrdft { @@ -37,6 +38,10 @@ class SpinLDA final : public Functional { ~SpinLDA() override = default; bool isSpin() const override { return true; } + bool isGGA() const override { return false; } + bool isMetaGGA() const override { return false; } + int numIn() const override { return 2; } + int numOut() const override { if (Factory::libxc) {return 3;} else {return xcfun_output_length(xcfun.get());} } private: mrcpp::FunctionTreeVector<3> rho_a; diff --git a/src/mrenv.cpp b/src/mrenv.cpp index 7cb93a412..e27417b3e 100644 --- a/src/mrenv.cpp +++ b/src/mrenv.cpp @@ -202,7 +202,6 @@ void mrenv::print_header() { print_utils::scalar(0, "Total cores ", (mrcpp::mpi::world_size - mrcpp::mpi::tot_bank_size) * mrcpp::omp::n_threads + mrcpp::mpi::tot_bank_size, "", 0, false); mrcpp::print::separator(0, ' '); mrcpp::print::separator(0, '-', 1); - printout(0, xcfun_splash()); mrcpp::print::environment(0); MRA->print(); } diff --git a/src/scf_solver/GroundStateSolver.cpp b/src/scf_solver/GroundStateSolver.cpp index a60985418..670ce8c53 100644 --- a/src/scf_solver/GroundStateSolver.cpp +++ b/src/scf_solver/GroundStateSolver.cpp @@ -194,6 +194,7 @@ void GroundStateSolver::printParameters(const std::string &calculation) const { mrcpp::print::separator(0, '~'); print_utils::text(0, "Calculation ", calculation); print_utils::text(0, "Method ", this->methodName); + print_utils::text(0, "XC library ", (this->libxc) ? "LibXC" : "XCFun"); print_utils::text(0, "Relativity ", this->relativityName); print_utils::text(0, "Environment ", this->environmentName); print_utils::text(0, "External fields ", this->externalFieldName); diff --git a/src/scf_solver/SCFSolver.h b/src/scf_solver/SCFSolver.h index 1c65b0f83..e54bcc236 100644 --- a/src/scf_solver/SCFSolver.h +++ b/src/scf_solver/SCFSolver.h @@ -57,6 +57,7 @@ class SCFSolver { void setHelmholtzPrec(double prec) { this->helmPrec = prec; } void setMaxIterations(int iter) { this->maxIter = iter; } void setMethodName(const std::string &name) { this->methodName = name; } + void setLibxc(bool libxc_) { this->libxc = libxc_; } void setRelativityName(const std::string &name) { this->relativityName = name; } void setEnvironmentName(const std::string &name) { this->environmentName = name; } void setExternalFieldName(const std::string &name) { this->externalFieldName = name; } @@ -73,6 +74,7 @@ class SCFSolver { std::string relativityName{"None"}; ///< Name of ZORA method std::string environmentName{"None"}; ///< Name for external environment std::string externalFieldName{"None"}; ///< Name for external fields + bool libxc; ///< Type of XC Library std::vector error; ///< Convergence orbital error std::vector property; ///< Convergence property error diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index df8fbc46a..e1aaea56a 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -30,10 +30,15 @@ add_subdirectory(h_el_field) add_subdirectory(h2_scf_hf) add_subdirectory(h2_pol_lda) add_subdirectory(h2_mag_lda) +add_subdirectory(h2_scf_svwn5_libxc) +add_subdirectory(h2_scf_b3lyp_libxc) add_subdirectory(h2o_energy_blyp) add_subdirectory(h2o_hirshfeld_lda) add_subdirectory(h2_pol_cube) +# TODO: Should rename test using xcfun to +"_xcfun" and add "xc_library = xcfun" add_subdirectory(li_scf_pbe0) +add_subdirectory(li_scf_b3lyp_libxc) +add_subdirectory(li_scf_svwn5_libxc) add_subdirectory(li_pol_lda) add_subdirectory(hf_grad_lda) add_subdirectory(hf_grad_blyp_surface_force) diff --git a/tests/h2_mag_lda/h2.inp b/tests/h2_mag_lda/h2.inp index e4389d90e..65b503dec 100644 --- a/tests/h2_mag_lda/h2.inp +++ b/tests/h2_mag_lda/h2.inp @@ -13,7 +13,7 @@ "restricted": false }, "DFT": { - "functionals": "LDA\n" + "functionals": "SVWN5\n" }, "Properties": { "magnetizability": true, diff --git a/tests/h2_pol_lda/h2.inp b/tests/h2_pol_lda/h2.inp index 0eab0fdb3..a59ae15a5 100644 --- a/tests/h2_pol_lda/h2.inp +++ b/tests/h2_pol_lda/h2.inp @@ -13,7 +13,7 @@ "restricted": false }, "DFT": { - "functionals": "LDA\n" + "functionals": "SVWN5\n" }, "Properties": { "polarizability": true diff --git a/tests/h2_scf_b3lyp_libxc/CMakeLists.txt b/tests/h2_scf_b3lyp_libxc/CMakeLists.txt new file mode 100644 index 000000000..c8c9dd34b --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/CMakeLists.txt @@ -0,0 +1,11 @@ +if(ENABLE_MPI) + set(_h2_scf_b3lyp_libxc_launcher "${MPIEXEC} ${MPIEXEC_NUMPROC_FLAG} 1") +endif() + +add_integration_test( + NAME "H2_SCF_B3LYP_Libxc" + LABELS "mrchem;H2_SCF_B3LYP_Libxc;libxc;dft;scf;energy" + COST 100 + LAUNCH_AGENT ${_h2_scf_b3lyp_libxc_launcher} + INITIAL_GUESS ${CMAKE_CURRENT_LIST_DIR}/initial_guess + ) diff --git a/tests/h2_scf_b3lyp_libxc/h2.inp b/tests/h2_scf_b3lyp_libxc/h2.inp new file mode 100644 index 000000000..6bf997644 --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/h2.inp @@ -0,0 +1,32 @@ +# vim:syntax=sh: + +world_prec = 1.0e-3 # Overall relative precision +world_size = 5 # Size of simulation box 2^n +world_unit = angstrom + +MPI { + numerically_exact = true # Guarantee identical results in MPI +} + +Molecule { +$coords +H 0.0000 0.0000 -0.3705 +H 0.0000 0.0000 0.3705 +$end +} + +WaveFunction { + method = b3lyp # Wave function method (HF or DFT) + restricted = true +} + +SCF { + kain = 3 # Length of KAIN iterative history + max_iter = 5 + orbital_thrs = 2.0e-2 + guess_type = GTO # Type of initial guess (none, gto, mw) +} + +DFT { + xc_library = libxc +} \ No newline at end of file diff --git a/tests/h2_scf_b3lyp_libxc/h2.out b/tests/h2_scf_b3lyp_libxc/h2.out new file mode 100644 index 000000000..eb6da44ad --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/h2.out @@ -0,0 +1,312 @@ + + +*************************************************************************** +*** *** +*** *** +*** __ __ ____ ____ _ *** +*** | \/ | _ \ / ___| |__ ___ _ __ ___ *** +*** | |\/| | |_) | | | '_ \ / _ \ '_ ` _ \ *** +*** | | | | _ <| |___| | | | __/ | | | | | *** +*** |_| |_|_| \_\\____|_| |_|\___|_| |_| |_| *** +*** *** +*** VERSION 1.2.0-alpha *** +*** *** +*** Git branch libxc_testbranch *** +*** Git commit hash 2bb9e2fa65c0876a30c8-dirty *** +*** Git commit author ylvao *** +*** Git commit date Thu Jan 29 12:29:18 2026 +0100 *** +*** *** +*** Contact: luca.frediani@uit.no *** +*** *** +*** Radovan Bast Magnar Bjorgve *** +*** Roberto Di Remigio Antoine Durdek *** +*** Luca Frediani Gabriel Gerez *** +*** Stig Rune Jensen Jonas Juselius *** +*** Rune Monstad Peter Wind *** +*** *** +*************************************************************************** + +--------------------------------------------------------------------------- + + MPI processes : (no bank) 1 + OpenMP threads : 1 + Total cores : 1 + +--------------------------------------------------------------------------- + + +--------------------------------------------------------------------------- + + MRCPP version : 1.6.0-alpha + Git branch : HEAD + Git commit hash : cec3cc96d578cdf77e8e + Git commit author : Niklas Göllmann + Git commit date : Tue Jan 6 16:29:30 2026 +0100 + + Linear algebra : EIGEN v3.4.0 + Parallelization : NONE + +--------------------------------------------------------------------------- + + + +=========================================================================== + MultiResolution Analysis +--------------------------------------------------------------------------- + polynomial order : 5 + polynomial type : Interpolating +--------------------------------------------------------------------------- + total boxes : 8 + boxes : [ 2 2 2 ] + unit lengths : [ 32.00000 32.00000 32.00000 ] + scaling factor : [ 1.00000 1.00000 1.00000 ] + lower bounds : [ -32.00000 -32.00000 -32.00000 ] + upper bounds : [ 32.00000 32.00000 32.00000 ] + total length : [ 64.00000 64.00000 64.00000 ] +=========================================================================== + + + +*************************************************************************** +*** *** +*** Initializing Molecule *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 1 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 H : 0.000000 0.000000 -0.700144 + 1 H : 0.000000 0.000000 0.700144 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +Name used in MapFunctionalName: B3LYP +Functional number: 0: b3lyp + + +*************************************************************************** +*** *** +*** XC Functional *** +*** *** +*************************************************************************** + +Using Libxc (version 7.0.0) to evaluate density functionals. Libxc is free +software. It is distributed under the Mozilla Public License, version 2.0. +For more information, please check the Libxc manual. You should cite + +S. Lehtola, C. Steigemann, M. J.T. Oliveira, and M. A.L. Marques., +SoftwareX 7, 1–5 (2018) DOI: 10.1016/j.softx.2017.11.002 + +when reporting the results of your calculation in a scientific article. + +Libxc functionals used in this calculation: + - hyb_gga_xc_b3lyp5 (ID 475): B3LYP with VWN functional 5 instead of RPA + * P. J. Stephens, F. J. Devlin, C. F. Chabalowski, and M. J. Frisch., J. Phys. Chem. 98, 11623 (1994) (DOI:10.1021/j100096a001) + +*************************************************************************** +*** *** +*** Computing Initial Guess Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial orbitals + Method : Project GTO molecular orbitals + Precision : 1.00000e-03 + Screening : 1.20000e+01 StdDev + Restricted : True + MO file : initial_guess/mrchem.mop +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Orbitals +--------------------------------------------------------------------------- + Alpha electrons : 1 + Beta electrons : 1 + Total electrons : 2 +--------------------------------------------------------------------------- + n Occ Spin : Norm +--------------------------------------------------------------------------- + 0 2 p : 9.999999041879e-01 +=========================================================================== + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial energy + Method : DFT (B3LYP) + XC Library : LibXC + Relativity : None + Environment : None + External fields : None + Precision : 1.00000e-03 + Density cutoff : 1.00000e-11 + Localization : Off +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +using libxc +=========================================================================== + Molecular Energy (initial) +--------------------------------------------------------------------------- + Kinetic energy : (au) 1.073486065232 + E-N energy : (au) -3.566631555860 + Coulomb energy : (au) 1.328363917978 + Exchange energy : (au) -0.132833922118 + X-C energy : (au) -0.558504025890 + N-N energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Electronic energy : (au) -1.856119520658 + Nuclear energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Total energy : (au) -1.141980234689e+00 + : (kcal/mol) -7.166034164599e+02 + : (kJ/mol) -2.998268694468e+03 + : (eV) -3.107486525140e+01 +=========================================================================== + + +=========================================================================== + Orbital Energies (initial) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 2 p : (au) -0.408773186422 +--------------------------------------------------------------------------- + Sum occupied : (au) -0.817546372845 +=========================================================================== + + + +*************************************************************************** +*** *** +*** Computing Ground State Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Optimize ground state orbitals + Method : DFT (B3LYP) + XC library : LibXC + Relativity : None + Environment : None + External fields : None + Checkpointing : Off + Max iterations : 5 + KAIN solver : 3 + Localization : Off + Diagonalization : First two iterations + Start precision : 1.00000e-03 + Final precision : 1.00000e-03 + Helmholtz precision : Dynamic + Energy threshold : Off + Orbital threshold : 2.00000e-02 +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Iter MO residual Total energy Update +--------------------------------------------------------------------------- + 0 1.000000e+00 -1.141980234689 -1.141980e+00 + 1 3.458102e-02 -1.148195369704 -6.215135e-03 + 2 4.033649e-03 -1.148449277513 -2.539078e-04 +--------------------------------------------------------------------------- + SCF converged in 2 iterations! +=========================================================================== + + + +*************************************************************************** +*** *** +*** Printing Molecular Properties *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 1 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 H : 0.000000 0.000000 -0.700144 + 1 H : 0.000000 0.000000 0.700144 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +=========================================================================== + Molecular Energy (final) +--------------------------------------------------------------------------- + Kinetic energy : (au) 1.062184984187 + E-N energy : (au) -3.541391438401 + Coulomb energy : (au) 1.294952081959 + Exchange energy : (au) -0.129492814104 + X-C energy : (au) -0.548841377124 + N-N energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Electronic energy : (au) -1.862588563482 + Nuclear energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Total energy : (au) -1.148449277513e+00 + : (kcal/mol) -7.206628021205e+02 + : (kJ/mol) -3.015253164072e+03 + : (eV) -3.125089687434e+01 +=========================================================================== + + +=========================================================================== + Orbital Energies (final) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 2 p : (au) -0.424856455858 +--------------------------------------------------------------------------- + Sum occupied : (au) -0.849712911716 +=========================================================================== + + +=========================================================================== + Dipole Moment (dip-1) +--------------------------------------------------------------------------- + r_O : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Electronic vector : -0.000000 0.000000 0.000004 + Magnitude : (au) 0.000004 + : (Debye) 0.000009 +--------------------------------------------------------------------------- + Nuclear vector : 0.000000 0.000000 -0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Total vector : 0.000000 0.000000 0.000004 + Magnitude : (au) 0.000004 + : (Debye) 0.000009 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** Exiting MRChem *** +*** *** +*** Wall time : 0h 0m 6s *** +*** *** +*************************************************************************** + + diff --git a/tests/h2_scf_b3lyp_libxc/initial_guess/mrchem.bas b/tests/h2_scf_b3lyp_libxc/initial_guess/mrchem.bas new file mode 100644 index 000000000..2b07b1622 --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/initial_guess/mrchem.bas @@ -0,0 +1,12 @@ +Gaussian basis cc-pVDZ + 1 + 1. 2 2 1 1 +H 0.0000000000 0.0000000000 -0.7000000000 +H 0.0000000000 0.0000000000 0.7000000000 + 4 2 + 13.0100000 0.01968500 0.00000000 + 1.9620000 0.13797700 0.00000000 + 0.4446000 0.47814800 0.00000000 + 0.1220000 0.50124000 1.00000000 + 1 1 + 0.7270000 1.00000000 diff --git a/tests/h2_scf_b3lyp_libxc/initial_guess/mrchem.mop b/tests/h2_scf_b3lyp_libxc/initial_guess/mrchem.mop new file mode 100644 index 000000000..fbc685ea2 --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/initial_guess/mrchem.mop @@ -0,0 +1,101 @@ + 10 + 0.679818844859533 + -0.160895416296538 + -0.000000000000000 + 0.000000000000000 + 0.017200409670392 + 0.679805451164590 + -0.160897270288160 + 0.000000000000000 + -0.000000000000000 + -0.017201747283775 + 0.394259467128243 + 1.587077397076630 + 0.000000000000000 + -0.000000000000000 + -0.022389357994487 + -0.394263910524846 + -1.587077979012975 + 0.000000000000000 + 0.000000000000000 + -0.022387397572700 + 1.187321376637329 + -1.316250637902104 + 0.000000000000000 + 0.000000000000000 + 0.022694692771659 + 1.187304649593490 + -1.316224276109744 + 0.000000000000000 + 0.000000000000000 + -0.022696339380126 + 1.376073054086166 + -2.492342381832811 + 0.000000000000000 + 0.000000000000000 + -0.358858997767801 + -1.376084250923773 + 2.492352716677140 + 0.000000000000000 + 0.000000000000000 + -0.358855133999848 + 0.000000000000000 + -0.000000000000001 + 0.154953780296779 + -0.558091716418119 + 0.000000000000000 + -0.000000000000000 + 0.000000000000001 + 0.154951569196481 + -0.558083752774001 + -0.000000000000000 + 0.000000000000002 + -0.000000000000002 + -0.558091716418118 + -0.154953780296779 + -0.000000000000000 + -0.000000000000001 + 0.000000000000001 + -0.558083752774003 + -0.154951569196481 + 0.000000000000000 + -0.768175736191248 + 0.618306842395530 + 0.000000000000002 + -0.000000000000000 + 0.726245917411092 + -0.768198542983742 + 0.618322374960124 + -0.000000000000002 + 0.000000000000001 + -0.726240941243852 + -0.000000000000001 + 0.000000000000001 + 0.197139927993266 + 0.970753591501604 + -0.000000000000000 + 0.000000000000001 + -0.000000000000000 + -0.197140889760886 + -0.970758327423854 + -0.000000000000001 + 0.000000000000000 + -0.000000000000001 + 0.970753591501605 + -0.197139927993266 + -0.000000000000002 + 0.000000000000003 + -0.000000000000002 + -0.970758327423854 + 0.197140889760886 + 0.000000000000001 + -4.523132495077299 + 2.174458822653225 + -0.000000000000001 + -0.000000000000000 + -2.033927630281827 + 4.523131231784129 + -2.174457955524351 + 0.000000000000000 + 0.000000000000000 + -2.033930080692429 diff --git a/tests/h2_scf_b3lyp_libxc/reference/h2.json b/tests/h2_scf_b3lyp_libxc/reference/h2.json new file mode 100644 index 000000000..9324d2438 --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/reference/h2.json @@ -0,0 +1,347 @@ +{ + "input": { + "constants": { + "N_a": 6.02214076e+23, + "angstrom2bohrs": 1.8897261246257702, + "boltzmann_constant": 1.380649e-23, + "dipmom_au2debye": 2.5417464739297717, + "e0": 8.8541878128e-12, + "electron_g_factor": -2.00231930436256, + "elementary_charge": 1.602176634e-19, + "fine_structure_constant": 0.0072973525693, + "hartree2ev": 27.211386245988, + "hartree2kcalmol": 627.5094740630558, + "hartree2kjmol": 2625.4996394798254, + "hartree2simagnetizability": 78.9451185, + "hartree2wavenumbers": 219474.6313632, + "light_speed": 137.035999084, + "meter2bohr": 18897261246.2577 + }, + "geom_opt": { + "init_step_size": -0.5, + "max_force_component": 0.005, + "max_history_length": 10, + "max_iter": 100, + "minimal_step_size": 0.01, + "run": false, + "subspace_tolerance": 0.001, + "use_previous_guess": false + }, + "molecule": { + "charge": 0, + "coords": [ + { + "atom": "h", + "r_rms": 2.6569547399e-05, + "xyz": [ + 0.0, + 0.0, + -0.7001435291738478 + ] + }, + { + "atom": "h", + "r_rms": 2.6569547399e-05, + "xyz": [ + 0.0, + 0.0, + 0.7001435291738478 + ] + } + ], + "multiplicity": 1 + }, + "mpi": { + "bank_per_node": -1, + "bank_size": -1, + "numerically_exact": true, + "omp_threads": -1, + "shared_memory_size": 10000, + "use_omp_num_threads": false + }, + "mra": { + "basis_order": 5, + "basis_type": "interpolating", + "boxes": [ + 2, + 2, + 2 + ], + "corner": [ + -1, + -1, + -1 + ], + "max_scale": 20, + "min_scale": -5 + }, + "printer": { + "file_name": "h2", + "print_constants": false, + "print_level": 0, + "print_mpi": false, + "print_prec": 6, + "print_width": 75 + }, + "rsp_calculations": {}, + "scf_calculation": { + "fock_operator": { + "coulomb_operator": { + "poisson_prec": 0.001, + "shared_memory": false + }, + "exchange_operator": { + "exchange_prec": -1.0, + "poisson_prec": 0.001 + }, + "kinetic_operator": { + "derivative": "abgv_55" + }, + "nuclear_operator": { + "nuclear_model": "point_like", + "proj_prec": 0.001, + "shared_memory": false, + "smooth_prec": 0.001 + }, + "xc_library": [ + "xcfun" + ], + "xc_operator": { + "shared_memory": false, + "xc_functional": { + "cutoff": 1e-11, + "functionals": [ + { + "coef": 1.0, + "name": "b3lyp" + } + ], + "spin": false + } + } + }, + "initial_guess": { + "cutoff": 1e-11, + "environment": "None", + "external_field": "None", + "file_CUBE_a": "cube_vectors/CUBE_a_vector.json", + "file_CUBE_b": "cube_vectors/CUBE_b_vector.json", + "file_CUBE_p": "cube_vectors/CUBE_p_vector.json", + "file_basis": "initial_guess/mrchem.bas", + "file_chk": "checkpoint/phi_scf", + "file_gto_a": "initial_guess/mrchem.moa", + "file_gto_b": "initial_guess/mrchem.mob", + "file_gto_p": "initial_guess/mrchem.mop", + "file_phi_a": "initial_guess/phi_a_scf", + "file_phi_b": "initial_guess/phi_b_scf", + "file_phi_p": "initial_guess/phi_p_scf", + "localize": false, + "method": "DFT", + "prec": 0.001, + "relativity": "None", + "restricted": true, + "rotate": true, + "screen": 12.0, + "type": "gto", + "xc_library": "xcfun", + "zeta": 0 + }, + "properties": { + "dipole_moment": { + "dip-1": { + "operator": "h_e_dip", + "precision": 0.001, + "r_O": [ + 0.0, + 0.0, + 0.0 + ] + } + } + }, + "scf_solver": { + "checkpoint": false, + "energy_thrs": -1.0, + "environment": "None", + "external_field": "None", + "file_chk": "checkpoint/phi_scf", + "final_prec": 0.001, + "helmholtz_prec": -1.0, + "kain": 3, + "localize": false, + "max_iter": 5, + "method": "DFT", + "orbital_thrs": 0.02, + "relativity": "None", + "rotation": 0, + "start_prec": 0.001, + "xc_library": "xcfun" + } + }, + "schema_name": "mrchem_input", + "schema_version": 1 + }, + "output": { + "properties": { + "center_of_mass": [ + 0.0, + 0.0, + 2.544077373851441e-17 + ], + "charge": 0, + "dipole_moment": { + "dip-1": { + "magnitude": 3.621717575627693e-06, + "r_O": [ + 0.0, + 0.0, + 0.0 + ], + "vector": [ + 0.0, + 0.0, + 3.6217175756276915e-06 + ], + "vector_el": [ + 0.0, + 0.0, + 3.621717584953565e-06 + ], + "vector_nuc": [ + 0.0, + 0.0, + 0.0 + ] + } + }, + "geometry": [ + { + "symbol": "H", + "xyz": [ + 0.0, + 0.0, + -0.7001435291738478 + ] + }, + { + "symbol": "H", + "xyz": [ + 0.0, + 0.0, + 0.7001435291738478 + ] + } + ], + "multiplicity": 1, + "orbital_energies": { + "energy": [ + -0.4248564815467927 + ], + "occupation": [ + 2.0 + ], + "spin": [ + "p" + ], + "sum_occupied": -0.8497129630935853 + }, + "scf_energy": { + "E_ee": 1.294952023161424, + "E_eext": 0.0, + "E_el": -1.8625885636832564, + "E_en": -3.541391351660119, + "E_kin": 1.0621849310130025, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.1484492777142954, + "E_x": -0.12949280822382017, + "E_xc": -0.5488413579737434, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + } + }, + "provenance": { + "creator": "MRChem", + "mpi_processes": 1, + "nthreads": 1, + "routine": "mrchem.x", + "total_cores": 1, + "version": "1.2.0-alpha" + }, + "rsp_calculations": null, + "scf_calculation": { + "initial_energy": { + "E_ee": 1.3283639179776123, + "E_eext": 0.0, + "E_el": -1.8561195208098742, + "E_en": -3.5666315558595945, + "E_kin": 1.0734860652323082, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.1419802348409132, + "E_x": -0.1328339221177218, + "E_xc": -0.5585040260424782, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "scf_solver": { + "converged": true, + "cycles": [ + { + "energy_terms": { + "E_ee": 1.2923397087862605, + "E_eext": 0.0, + "E_el": -1.8623346560358298, + "E_en": -3.5288499596074367, + "E_kin": 1.050814208659165, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.148195370066869, + "E_x": -0.1292315830207368, + "E_xc": -0.5474070308530816, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -1.148195370066869, + "energy_update": 0.006215135225955848, + "mo_residual": 0.03458114843901992, + "wall_time": 5.117739776 + }, + { + "energy_terms": { + "E_ee": 1.294952023161424, + "E_eext": 0.0, + "E_el": -1.8625885636832564, + "E_en": -3.541391351660119, + "E_kin": 1.0621849310130025, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.1484492777142954, + "E_x": -0.12949280822382017, + "E_xc": -0.5488413579737434, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -1.1484492777142954, + "energy_update": 0.00025390764742638794, + "mo_residual": 0.004033663714494966, + "wall_time": 3.594082832 + } + ], + "wall_time": 8.713920775 + }, + "success": true + }, + "schema_name": "mrchem_output", + "schema_version": 1, + "success": true + } +} diff --git a/tests/h2_scf_b3lyp_libxc/reference/h2.out b/tests/h2_scf_b3lyp_libxc/reference/h2.out new file mode 100644 index 000000000..1affaefaf --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/reference/h2.out @@ -0,0 +1,310 @@ + + +*************************************************************************** +*** *** +*** *** +*** __ __ ____ ____ _ *** +*** | \/ | _ \ / ___| |__ ___ _ __ ___ *** +*** | |\/| | |_) | | | '_ \ / _ \ '_ ` _ \ *** +*** | | | | _ <| |___| | | | __/ | | | | | *** +*** |_| |_|_| \_\\____|_| |_|\___|_| |_| |_| *** +*** *** +*** VERSION 1.2.0-alpha *** +*** *** +*** Git branch libxc_testbranch *** +*** Git commit hash 2bb9e2fa65c0876a30c8 *** +*** Git commit author ylvao *** +*** Git commit date Thu Jan 29 12:29:18 2026 +0100 *** +*** *** +*** Contact: luca.frediani@uit.no *** +*** *** +*** Radovan Bast Magnar Bjorgve *** +*** Roberto Di Remigio Antoine Durdek *** +*** Luca Frediani Gabriel Gerez *** +*** Stig Rune Jensen Jonas Juselius *** +*** Rune Monstad Peter Wind *** +*** *** +*************************************************************************** + +--------------------------------------------------------------------------- + + MPI processes : (no bank) 1 + OpenMP threads : 1 + Total cores : 1 + +--------------------------------------------------------------------------- + + +--------------------------------------------------------------------------- + + MRCPP version : 1.6.0-alpha + Git branch : HEAD + Git commit hash : cec3cc96d578cdf77e8e + Git commit author : Niklas Göllmann + Git commit date : Tue Jan 6 16:29:30 2026 +0100 + + Linear algebra : EIGEN v3.4.0 + Parallelization : NONE + +--------------------------------------------------------------------------- + + + +=========================================================================== + MultiResolution Analysis +--------------------------------------------------------------------------- + polynomial order : 5 + polynomial type : Interpolating +--------------------------------------------------------------------------- + total boxes : 8 + boxes : [ 2 2 2 ] + unit lengths : [ 32.00000 32.00000 32.00000 ] + scaling factor : [ 1.00000 1.00000 1.00000 ] + lower bounds : [ -32.00000 -32.00000 -32.00000 ] + upper bounds : [ 32.00000 32.00000 32.00000 ] + total length : [ 64.00000 64.00000 64.00000 ] +=========================================================================== + + + +*************************************************************************** +*** *** +*** Initializing Molecule *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 1 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 H : 0.000000 0.000000 -0.700144 + 1 H : 0.000000 0.000000 0.700144 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** XC Functional *** +*** *** +*************************************************************************** + +XCFun DFT library Copyright 2009-2020 Ulf Ekstrom and contributors. +See http://dftlibs.org/xcfun/ for more information. + +This is free software; see the source code for copying conditions. +There is ABSOLUTELY NO WARRANTY; not even for MERCHANTABILITY or +FITNESS FOR A PARTICULAR PURPOSE. For details see the documentation. +Scientific users of this library should cite +U. Ekstrom, L. Visscher, R. Bast, A. J. Thorvaldsen and K. Ruud; +J.Chem.Theor.Comp. 2010, DOI: 10.1021/ct100117s + +--------------------------------------------------------------------------- + + +*************************************************************************** +*** *** +*** Computing Initial Guess Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial orbitals + Method : Project GTO molecular orbitals + Precision : 1.00000e-03 + Screening : 1.20000e+01 StdDev + Restricted : True + MO file : initial_guess/mrchem.mop +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Orbitals +--------------------------------------------------------------------------- + Alpha electrons : 1 + Beta electrons : 1 + Total electrons : 2 +--------------------------------------------------------------------------- + n Occ Spin : Norm +--------------------------------------------------------------------------- + 0 2 p : 9.999999041879e-01 +=========================================================================== + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial energy + Method : DFT + XC Library : XCFun + Relativity : None + Environment : None + External fields : None + Precision : 1.00000e-03 + Density cutoff : 1.00000e-11 + Localization : Off +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +using xcfun +=========================================================================== + Molecular Energy (initial) +--------------------------------------------------------------------------- + Kinetic energy : (au) 1.073486065232 + E-N energy : (au) -3.566631555860 + Coulomb energy : (au) 1.328363917978 + Exchange energy : (au) -0.132833922118 + X-C energy : (au) -0.558504026042 + N-N energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Electronic energy : (au) -1.856119520810 + Nuclear energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Total energy : (au) -1.141980234841e+00 + : (kcal/mol) -7.166034165554e+02 + : (kJ/mol) -2.998268694868e+03 + : (eV) -3.107486525554e+01 +=========================================================================== + + +=========================================================================== + Orbital Energies (initial) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 2 p : (au) -0.408773186621 +--------------------------------------------------------------------------- + Sum occupied : (au) -0.817546373242 +=========================================================================== + + + +*************************************************************************** +*** *** +*** Computing Ground State Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Optimize ground state orbitals + Method : DFT + XC library : XCFun + Relativity : None + Environment : None + External fields : None + Checkpointing : Off + Max iterations : 5 + KAIN solver : 3 + Localization : Off + Diagonalization : First two iterations + Start precision : 1.00000e-03 + Final precision : 1.00000e-03 + Helmholtz precision : Dynamic + Energy threshold : Off + Orbital threshold : 2.00000e-02 +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Iter MO residual Total energy Update +--------------------------------------------------------------------------- + 0 1.000000e+00 -1.141980234841 -1.141980e+00 + 1 3.458115e-02 -1.148195370067 -6.215135e-03 + 2 4.033664e-03 -1.148449277714 -2.539076e-04 +--------------------------------------------------------------------------- + SCF converged in 2 iterations! +=========================================================================== + + + +*************************************************************************** +*** *** +*** Printing Molecular Properties *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 1 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 H : 0.000000 0.000000 -0.700144 + 1 H : 0.000000 0.000000 0.700144 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +=========================================================================== + Molecular Energy (final) +--------------------------------------------------------------------------- + Kinetic energy : (au) 1.062184931013 + E-N energy : (au) -3.541391351660 + Coulomb energy : (au) 1.294952023161 + Exchange energy : (au) -0.129492808224 + X-C energy : (au) -0.548841357974 + N-N energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Electronic energy : (au) -1.862588563683 + Nuclear energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Total energy : (au) -1.148449277714e+00 + : (kcal/mol) -7.206628022466e+02 + : (kJ/mol) -3.015253164600e+03 + : (eV) -3.125089687981e+01 +=========================================================================== + + +=========================================================================== + Orbital Energies (final) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 2 p : (au) -0.424856481547 +--------------------------------------------------------------------------- + Sum occupied : (au) -0.849712963094 +=========================================================================== + + +=========================================================================== + Dipole Moment (dip-1) +--------------------------------------------------------------------------- + r_O : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Electronic vector : -0.000000 -0.000000 0.000004 + Magnitude : (au) 0.000004 + : (Debye) 0.000009 +--------------------------------------------------------------------------- + Nuclear vector : 0.000000 0.000000 -0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Total vector : 0.000000 0.000000 0.000004 + Magnitude : (au) 0.000004 + : (Debye) 0.000009 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** Exiting MRChem *** +*** *** +*** Wall time : 0h 0m 10s *** +*** *** +*************************************************************************** + + diff --git a/tests/h2_scf_b3lyp_libxc/test b/tests/h2_scf_b3lyp_libxc/test new file mode 100755 index 000000000..01254e21c --- /dev/null +++ b/tests/h2_scf_b3lyp_libxc/test @@ -0,0 +1,26 @@ +#!/usr/bin/env python3 + +import sys +from pathlib import Path + +sys.path.append(str(Path(__file__).resolve().parents[1])) + +from tester import * # isort:skip + +options = script_cli() + +filters = { + SUM_OCCUPIED: rel_tolerance(1.0e-7), + E_KIN: rel_tolerance(1.0e-7), + E_EN: rel_tolerance(1.0e-7), + E_EE: rel_tolerance(1.0e-7), + E_X: rel_tolerance(1.0e-7), + E_XC: rel_tolerance(1.0e-7), + E_EEXT: rel_tolerance(1.0e-7), + E_NEXT: rel_tolerance(1.0e-7), + E_EL: rel_tolerance(1.0e-7), +} + +ierr = run(options, input_file="h2", filters=filters) + +sys.exit(ierr) diff --git a/tests/h2_scf_lda_libxc/h2.inp b/tests/h2_scf_lda_libxc/h2.inp new file mode 100644 index 000000000..4cc3bca1b --- /dev/null +++ b/tests/h2_scf_lda_libxc/h2.inp @@ -0,0 +1,32 @@ +# vim:syntax=sh: + +world_prec = 1.0e-3 # Overall relative precision +world_size = 5 # Size of simulation box 2^n +world_unit = angstrom + +MPI { + numerically_exact = true # Guarantee identical results in MPI +} + +Molecule { +$coords +H 0.0000 0.0000 -0.3705 +H 0.0000 0.0000 0.3705 +$end +} + +WaveFunction { + method = SVWN5 # Wave function method (HF or DFT) + restricted = true +} + +SCF { + kain = 3 # Length of KAIN iterative history + max_iter = 5 + orbital_thrs = 2.0e-2 + guess_type = GTO # Type of initial guess (none, gto, mw) +} + +DFT { + xc_library = libxc +} \ No newline at end of file diff --git a/tests/h2_scf_svwn5_libxc/CMakeLists.txt b/tests/h2_scf_svwn5_libxc/CMakeLists.txt new file mode 100644 index 000000000..4f74ae912 --- /dev/null +++ b/tests/h2_scf_svwn5_libxc/CMakeLists.txt @@ -0,0 +1,11 @@ +if(ENABLE_MPI) + set(_h2_scf_svwn5_libxc "${MPIEXEC} ${MPIEXEC_NUMPROC_FLAG} 1") +endif() + +add_integration_test( + NAME "H2_SCF_SVWN5_Libxc" + LABELS "mrchem;h2_scf_svwn5_libxc;libxc;H2_SCF_SVWN5_Libxc;dft;scf;energy" + COST 100 + LAUNCH_AGENT ${_h2_scf_svwn5_libxc_launcher} + INITIAL_GUESS ${CMAKE_CURRENT_LIST_DIR}/initial_guess + ) diff --git a/tests/h2_scf_svwn5_libxc/h2.inp b/tests/h2_scf_svwn5_libxc/h2.inp new file mode 100644 index 000000000..a67075e2f --- /dev/null +++ b/tests/h2_scf_svwn5_libxc/h2.inp @@ -0,0 +1,32 @@ +# vim:syntax=sh: + +world_prec = 1.0e-3 # Overall relative precision +world_size = 5 # Size of simulation box 2^n +world_unit = angstrom + +MPI { + numerically_exact = true # Guarantee identical results in MPI +} + +Molecule { +$coords +H 0.0000 0.0000 -0.3705 +H 0.0000 0.0000 0.3705 +$end +} + +WaveFunction { + method = svwn5 # Wave function method (HF or DFT) + restricted = true +} + +SCF { + kain = 3 # Length of KAIN iterative history + max_iter = 5 + orbital_thrs = 2.0e-2 + guess_type = GTO # Type of initial guess (none, gto, mw) +} + +DFT { + xc_library = libxc +} \ No newline at end of file diff --git a/tests/h2_scf_svwn5_libxc/initial_guess/mrchem.bas b/tests/h2_scf_svwn5_libxc/initial_guess/mrchem.bas new file mode 100644 index 000000000..2b07b1622 --- /dev/null +++ b/tests/h2_scf_svwn5_libxc/initial_guess/mrchem.bas @@ -0,0 +1,12 @@ +Gaussian basis cc-pVDZ + 1 + 1. 2 2 1 1 +H 0.0000000000 0.0000000000 -0.7000000000 +H 0.0000000000 0.0000000000 0.7000000000 + 4 2 + 13.0100000 0.01968500 0.00000000 + 1.9620000 0.13797700 0.00000000 + 0.4446000 0.47814800 0.00000000 + 0.1220000 0.50124000 1.00000000 + 1 1 + 0.7270000 1.00000000 diff --git a/tests/h2_scf_svwn5_libxc/initial_guess/mrchem.mop b/tests/h2_scf_svwn5_libxc/initial_guess/mrchem.mop new file mode 100644 index 000000000..fbc685ea2 --- /dev/null +++ b/tests/h2_scf_svwn5_libxc/initial_guess/mrchem.mop @@ -0,0 +1,101 @@ + 10 + 0.679818844859533 + -0.160895416296538 + -0.000000000000000 + 0.000000000000000 + 0.017200409670392 + 0.679805451164590 + -0.160897270288160 + 0.000000000000000 + -0.000000000000000 + -0.017201747283775 + 0.394259467128243 + 1.587077397076630 + 0.000000000000000 + -0.000000000000000 + -0.022389357994487 + -0.394263910524846 + -1.587077979012975 + 0.000000000000000 + 0.000000000000000 + -0.022387397572700 + 1.187321376637329 + -1.316250637902104 + 0.000000000000000 + 0.000000000000000 + 0.022694692771659 + 1.187304649593490 + -1.316224276109744 + 0.000000000000000 + 0.000000000000000 + -0.022696339380126 + 1.376073054086166 + -2.492342381832811 + 0.000000000000000 + 0.000000000000000 + -0.358858997767801 + -1.376084250923773 + 2.492352716677140 + 0.000000000000000 + 0.000000000000000 + -0.358855133999848 + 0.000000000000000 + -0.000000000000001 + 0.154953780296779 + -0.558091716418119 + 0.000000000000000 + -0.000000000000000 + 0.000000000000001 + 0.154951569196481 + -0.558083752774001 + -0.000000000000000 + 0.000000000000002 + -0.000000000000002 + -0.558091716418118 + -0.154953780296779 + -0.000000000000000 + -0.000000000000001 + 0.000000000000001 + -0.558083752774003 + -0.154951569196481 + 0.000000000000000 + -0.768175736191248 + 0.618306842395530 + 0.000000000000002 + -0.000000000000000 + 0.726245917411092 + -0.768198542983742 + 0.618322374960124 + -0.000000000000002 + 0.000000000000001 + -0.726240941243852 + -0.000000000000001 + 0.000000000000001 + 0.197139927993266 + 0.970753591501604 + -0.000000000000000 + 0.000000000000001 + -0.000000000000000 + -0.197140889760886 + -0.970758327423854 + -0.000000000000001 + 0.000000000000000 + -0.000000000000001 + 0.970753591501605 + -0.197139927993266 + -0.000000000000002 + 0.000000000000003 + -0.000000000000002 + -0.970758327423854 + 0.197140889760886 + 0.000000000000001 + -4.523132495077299 + 2.174458822653225 + -0.000000000000001 + -0.000000000000000 + -2.033927630281827 + 4.523131231784129 + -2.174457955524351 + 0.000000000000000 + 0.000000000000000 + -2.033930080692429 diff --git a/tests/h2_scf_svwn5_libxc/reference/h2.json b/tests/h2_scf_svwn5_libxc/reference/h2.json new file mode 100644 index 000000000..2b33053b7 --- /dev/null +++ b/tests/h2_scf_svwn5_libxc/reference/h2.json @@ -0,0 +1,347 @@ +{ + "input": { + "constants": { + "N_a": 6.02214076e+23, + "angstrom2bohrs": 1.8897261246257702, + "boltzmann_constant": 1.380649e-23, + "dipmom_au2debye": 2.5417464739297717, + "e0": 8.8541878128e-12, + "electron_g_factor": -2.00231930436256, + "elementary_charge": 1.602176634e-19, + "fine_structure_constant": 0.0072973525693, + "hartree2ev": 27.211386245988, + "hartree2kcalmol": 627.5094740630558, + "hartree2kjmol": 2625.4996394798254, + "hartree2simagnetizability": 78.9451185, + "hartree2wavenumbers": 219474.6313632, + "light_speed": 137.035999084, + "meter2bohr": 18897261246.2577 + }, + "geom_opt": { + "init_step_size": -0.5, + "max_force_component": 0.005, + "max_history_length": 10, + "max_iter": 100, + "minimal_step_size": 0.01, + "run": false, + "subspace_tolerance": 0.001, + "use_previous_guess": false + }, + "molecule": { + "charge": 0, + "coords": [ + { + "atom": "h", + "r_rms": 2.6569547399e-05, + "xyz": [ + 0.0, + 0.0, + -0.7001435291738478 + ] + }, + { + "atom": "h", + "r_rms": 2.6569547399e-05, + "xyz": [ + 0.0, + 0.0, + 0.7001435291738478 + ] + } + ], + "multiplicity": 1 + }, + "mpi": { + "bank_per_node": -1, + "bank_size": -1, + "numerically_exact": true, + "omp_threads": -1, + "shared_memory_size": 10000, + "use_omp_num_threads": false + }, + "mra": { + "basis_order": 5, + "basis_type": "interpolating", + "boxes": [ + 2, + 2, + 2 + ], + "corner": [ + -1, + -1, + -1 + ], + "max_scale": 20, + "min_scale": -5 + }, + "printer": { + "file_name": "h2", + "print_constants": false, + "print_level": 0, + "print_mpi": false, + "print_prec": 6, + "print_width": 75 + }, + "rsp_calculations": {}, + "scf_calculation": { + "fock_operator": { + "coulomb_operator": { + "poisson_prec": 0.001, + "shared_memory": false + }, + "exchange_operator": { + "exchange_prec": -1.0, + "poisson_prec": 0.001 + }, + "kinetic_operator": { + "derivative": "abgv_55" + }, + "nuclear_operator": { + "nuclear_model": "point_like", + "proj_prec": 0.001, + "shared_memory": false, + "smooth_prec": 0.001 + }, + "xc_library": [ + "xcfun" + ], + "xc_operator": { + "shared_memory": false, + "xc_functional": { + "cutoff": 1e-11, + "functionals": [ + { + "coef": 1.0, + "name": "svwn5" + } + ], + "spin": false + } + } + }, + "initial_guess": { + "cutoff": 1e-11, + "environment": "None", + "external_field": "None", + "file_CUBE_a": "cube_vectors/CUBE_a_vector.json", + "file_CUBE_b": "cube_vectors/CUBE_b_vector.json", + "file_CUBE_p": "cube_vectors/CUBE_p_vector.json", + "file_basis": "initial_guess/mrchem.bas", + "file_chk": "checkpoint/phi_scf", + "file_gto_a": "initial_guess/mrchem.moa", + "file_gto_b": "initial_guess/mrchem.mob", + "file_gto_p": "initial_guess/mrchem.mop", + "file_phi_a": "initial_guess/phi_a_scf", + "file_phi_b": "initial_guess/phi_b_scf", + "file_phi_p": "initial_guess/phi_p_scf", + "localize": false, + "method": "DFT (SVWN5)", + "prec": 0.001, + "relativity": "None", + "restricted": true, + "rotate": true, + "screen": 12.0, + "type": "gto", + "xc_library": "xcfun", + "zeta": 0 + }, + "properties": { + "dipole_moment": { + "dip-1": { + "operator": "h_e_dip", + "precision": 0.001, + "r_O": [ + 0.0, + 0.0, + 0.0 + ] + } + } + }, + "scf_solver": { + "checkpoint": false, + "energy_thrs": -1.0, + "environment": "None", + "external_field": "None", + "file_chk": "checkpoint/phi_scf", + "final_prec": 0.001, + "helmholtz_prec": -1.0, + "kain": 3, + "localize": false, + "max_iter": 5, + "method": "DFT (SVWN5)", + "orbital_thrs": 0.02, + "relativity": "None", + "rotation": 0, + "start_prec": 0.001, + "xc_library": "xcfun" + } + }, + "schema_name": "mrchem_input", + "schema_version": 1 + }, + "output": { + "properties": { + "center_of_mass": [ + 0.0, + 0.0, + 2.544077373851441e-17 + ], + "charge": 0, + "dipole_moment": { + "dip-1": { + "magnitude": 2.9397596776633593e-06, + "r_O": [ + 0.0, + 0.0, + 0.0 + ], + "vector": [ + 0.0, + 0.0, + 2.939759677663358e-06 + ], + "vector_el": [ + 0.0, + 0.0, + 2.9397596869892314e-06 + ], + "vector_nuc": [ + 0.0, + 0.0, + 0.0 + ] + } + }, + "geometry": [ + { + "symbol": "H", + "xyz": [ + 0.0, + 0.0, + -0.7001435291738478 + ] + }, + { + "symbol": "H", + "xyz": [ + 0.0, + 0.0, + 0.7001435291738478 + ] + } + ], + "multiplicity": 1, + "orbital_energies": { + "energy": [ + -0.37321364059845763 + ], + "occupation": [ + 2.0 + ], + "spin": [ + "p" + ], + "sum_occupied": -0.7464272811969153 + }, + "scf_energy": { + "E_ee": 1.2731202320872728, + "E_eext": 0.0, + "E_el": -1.8260610902621828, + "E_en": -3.4931898634358642, + "E_kin": 1.0262233331403021, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.111921804293222, + "E_x": 0.0, + "E_xc": -0.6322147920538933, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + } + }, + "provenance": { + "creator": "MRChem", + "mpi_processes": 1, + "nthreads": 64, + "routine": "mrchem.x", + "total_cores": 64, + "version": "1.2.0-alpha" + }, + "rsp_calculations": null, + "scf_calculation": { + "initial_energy": { + "E_ee": 1.3283639179776123, + "E_eext": 0.0, + "E_el": -1.8186520281572673, + "E_en": -3.5666315558595945, + "E_kin": 1.0734860652323082, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.1045127421883065, + "E_x": 0.0, + "E_xc": -0.6538704555075932, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "scf_solver": { + "converged": true, + "cycles": [ + { + "energy_terms": { + "E_ee": 1.2666405198177648, + "E_eext": 0.0, + "E_el": -1.8258862964628415, + "E_en": -3.4756075194966325, + "E_kin": 1.0123072761988723, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.1117470104938807, + "E_x": 0.0, + "E_xc": -0.6292265729828465, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -1.1117470104938807, + "energy_update": 0.007234268305574165, + "mo_residual": 0.04926458544567805, + "wall_time": 1.957428272 + }, + { + "energy_terms": { + "E_ee": 1.2731202320872728, + "E_eext": 0.0, + "E_el": -1.8260610902621828, + "E_en": -3.4931898634358642, + "E_kin": 1.0262233331403021, + "E_next": 0.0, + "E_nn": 0.7141392859689609, + "E_nuc": 0.7141392859689609, + "E_tot": -1.111921804293222, + "E_x": 0.0, + "E_xc": -0.6322147920538933, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -1.111921804293222, + "energy_update": 0.00017479379934126626, + "mo_residual": 0.0062232191738004535, + "wall_time": 2.759983648 + } + ], + "wall_time": 4.719531115 + }, + "success": true + }, + "schema_name": "mrchem_output", + "schema_version": 1, + "success": true + } +} diff --git a/tests/h2_scf_svwn5_libxc/reference/h2.out b/tests/h2_scf_svwn5_libxc/reference/h2.out new file mode 100644 index 000000000..c8979d857 --- /dev/null +++ b/tests/h2_scf_svwn5_libxc/reference/h2.out @@ -0,0 +1,309 @@ + + +*************************************************************************** +*** *** +*** *** +*** __ __ ____ ____ _ *** +*** | \/ | _ \ / ___| |__ ___ _ __ ___ *** +*** | |\/| | |_) | | | '_ \ / _ \ '_ ` _ \ *** +*** | | | | _ <| |___| | | | __/ | | | | | *** +*** |_| |_|_| \_\\____|_| |_|\___|_| |_| |_| *** +*** *** +*** VERSION 1.2.0-alpha *** +*** *** +*** Git branch libxc_testbranch *** +*** Git commit hash cdd0549cf55900c1b4d7 *** +*** Git commit author ylvao *** +*** Git commit date Mon Feb 2 15:24:40 2026 +0100 *** +*** *** +*** Contact: luca.frediani@uit.no *** +*** *** +*** Radovan Bast Magnar Bjorgve *** +*** Roberto Di Remigio Antoine Durdek *** +*** Luca Frediani Gabriel Gerez *** +*** Stig Rune Jensen Jonas Juselius *** +*** Rune Monstad Peter Wind *** +*** *** +*************************************************************************** + +--------------------------------------------------------------------------- + + MPI processes : (no bank) 1 + OpenMP threads : 64 + Total cores : 64 + +--------------------------------------------------------------------------- + + +--------------------------------------------------------------------------- + + MRCPP version : 1.6.0-alpha + Git branch : HEAD + Git commit hash : cec3cc96d578cdf77e8e + Git commit author : Niklas Göllmann + Git commit date : Tue Jan 6 16:29:30 2026 +0100 + + Linear algebra : EIGEN v3.4.0 + Parallelization : OpenMP (64 threads) + +--------------------------------------------------------------------------- + + + +=========================================================================== + MultiResolution Analysis +--------------------------------------------------------------------------- + polynomial order : 5 + polynomial type : Interpolating +--------------------------------------------------------------------------- + total boxes : 8 + boxes : [ 2 2 2 ] + unit lengths : [ 32.00000 32.00000 32.00000 ] + scaling factor : [ 1.00000 1.00000 1.00000 ] + lower bounds : [ -32.00000 -32.00000 -32.00000 ] + upper bounds : [ 32.00000 32.00000 32.00000 ] + total length : [ 64.00000 64.00000 64.00000 ] +=========================================================================== + + + +*************************************************************************** +*** *** +*** Initializing Molecule *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 1 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 H : 0.000000 0.000000 -0.700144 + 1 H : 0.000000 0.000000 0.700144 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** XC Functional *** +*** *** +*************************************************************************** + +XCFun DFT library Copyright 2009-2020 Ulf Ekstrom and contributors. +See http://dftlibs.org/xcfun/ for more information. + +This is free software; see the source code for copying conditions. +There is ABSOLUTELY NO WARRANTY; not even for MERCHANTABILITY or +FITNESS FOR A PARTICULAR PURPOSE. For details see the documentation. +Scientific users of this library should cite +U. Ekstrom, L. Visscher, R. Bast, A. J. Thorvaldsen and K. Ruud; +J.Chem.Theor.Comp. 2010, DOI: 10.1021/ct100117s + +--------------------------------------------------------------------------- + + +*************************************************************************** +*** *** +*** Computing Initial Guess Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial orbitals + Method : Project GTO molecular orbitals + Precision : 1.00000e-03 + Screening : 1.20000e+01 StdDev + Restricted : True + MO file : initial_guess/mrchem.mop +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Orbitals +--------------------------------------------------------------------------- + Alpha electrons : 1 + Beta electrons : 1 + Total electrons : 2 +--------------------------------------------------------------------------- + n Occ Spin : Norm +--------------------------------------------------------------------------- + 0 2 p : 9.999999041879e-01 +=========================================================================== + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial energy + Method : DFT (SVWN5) + XC Library : XCFun + Relativity : None + Environment : None + External fields : None + Precision : 1.00000e-03 + Density cutoff : 1.00000e-11 + Localization : Off +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Energy (initial) +--------------------------------------------------------------------------- + Kinetic energy : (au) 1.073486065232 + E-N energy : (au) -3.566631555860 + Coulomb energy : (au) 1.328363917978 + Exchange energy : (au) 0.000000000000 + X-C energy : (au) -0.653870455508 + N-N energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Electronic energy : (au) -1.818652028157 + Nuclear energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Total energy : (au) -1.104512742188e+00 + : (kcal/mol) -6.930922099465e+02 + : (kJ/mol) -2.899897806416e+03 + : (eV) -3.005532284130e+01 +=========================================================================== + + +=========================================================================== + Orbital Energies (initial) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 2 p : (au) -0.345283968279 +--------------------------------------------------------------------------- + Sum occupied : (au) -0.690567936558 +=========================================================================== + + + +*************************************************************************** +*** *** +*** Computing Ground State Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Optimize ground state orbitals + Method : DFT (SVWN5) + XC library : XCFun + Relativity : None + Environment : None + External fields : None + Checkpointing : Off + Max iterations : 5 + KAIN solver : 3 + Localization : Off + Diagonalization : First two iterations + Start precision : 1.00000e-03 + Final precision : 1.00000e-03 + Helmholtz precision : Dynamic + Energy threshold : Off + Orbital threshold : 2.00000e-02 +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Iter MO residual Total energy Update +--------------------------------------------------------------------------- + 0 1.000000e+00 -1.104512742188 -1.104513e+00 + 1 4.926459e-02 -1.111747010494 -7.234268e-03 + 2 6.223219e-03 -1.111921804293 -1.747938e-04 +--------------------------------------------------------------------------- + SCF converged in 2 iterations! +=========================================================================== + + + +*************************************************************************** +*** *** +*** Printing Molecular Properties *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 1 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 H : 0.000000 0.000000 -0.700144 + 1 H : 0.000000 0.000000 0.700144 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +=========================================================================== + Molecular Energy (final) +--------------------------------------------------------------------------- + Kinetic energy : (au) 1.026223333140 + E-N energy : (au) -3.493189863436 + Coulomb energy : (au) 1.273120232087 + Exchange energy : (au) 0.000000000000 + X-C energy : (au) -0.632214792054 + N-N energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Electronic energy : (au) -1.826061090262 + Nuclear energy : (au) 0.714139285969 +--------------------------------------------------------------------------- + Total energy : (au) -1.111921804293e+00 + : (kcal/mol) -6.977414666113e+02 + : (kJ/mol) -2.919350296302e+03 + : (eV) -3.025693369196e+01 +=========================================================================== + + +=========================================================================== + Orbital Energies (final) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 2 p : (au) -0.373213640598 +--------------------------------------------------------------------------- + Sum occupied : (au) -0.746427281197 +=========================================================================== + + +=========================================================================== + Dipole Moment (dip-1) +--------------------------------------------------------------------------- + r_O : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Electronic vector : -0.000000 -0.000000 0.000003 + Magnitude : (au) 0.000003 + : (Debye) 0.000007 +--------------------------------------------------------------------------- + Nuclear vector : 0.000000 0.000000 -0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Total vector : 0.000000 0.000000 0.000003 + Magnitude : (au) 0.000003 + : (Debye) 0.000007 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** Exiting MRChem *** +*** *** +*** Wall time : 0h 0m 6s *** +*** *** +*************************************************************************** + + diff --git a/tests/h2_scf_svwn5_libxc/test b/tests/h2_scf_svwn5_libxc/test new file mode 100755 index 000000000..01254e21c --- /dev/null +++ b/tests/h2_scf_svwn5_libxc/test @@ -0,0 +1,26 @@ +#!/usr/bin/env python3 + +import sys +from pathlib import Path + +sys.path.append(str(Path(__file__).resolve().parents[1])) + +from tester import * # isort:skip + +options = script_cli() + +filters = { + SUM_OCCUPIED: rel_tolerance(1.0e-7), + E_KIN: rel_tolerance(1.0e-7), + E_EN: rel_tolerance(1.0e-7), + E_EE: rel_tolerance(1.0e-7), + E_X: rel_tolerance(1.0e-7), + E_XC: rel_tolerance(1.0e-7), + E_EEXT: rel_tolerance(1.0e-7), + E_NEXT: rel_tolerance(1.0e-7), + E_EL: rel_tolerance(1.0e-7), +} + +ierr = run(options, input_file="h2", filters=filters) + +sys.exit(ierr) diff --git a/tests/he_azora_scf_lda/he.inp b/tests/he_azora_scf_lda/he.inp index 3e1d0da66..f97688083 100644 --- a/tests/he_azora_scf_lda/he.inp +++ b/tests/he_azora_scf_lda/he.inp @@ -1,7 +1,7 @@ world_prec = 1.0e-4 WaveFunction { - method = LDA + method = SVWN5 restricted = true relativity = azora } diff --git a/tests/hf_grad_lda/hf.inp b/tests/hf_grad_lda/hf.inp index 43ebe5f04..ebee9a96d 100644 --- a/tests/hf_grad_lda/hf.inp +++ b/tests/hf_grad_lda/hf.inp @@ -19,7 +19,7 @@ Forces { } WaveFunction { - method = LDA # Wave function method (HF or DFT) + method = SVWN5 # Wave function method (HF or DFT) } Properties { diff --git a/tests/li_pol_lda/li.inp b/tests/li_pol_lda/li.inp index c09eef022..bd3ccfc69 100644 --- a/tests/li_pol_lda/li.inp +++ b/tests/li_pol_lda/li.inp @@ -15,7 +15,7 @@ $end } WaveFunction { - method = LDA # Wave function method (HF or DFT) + method = SVWN5 # Wave function method (HF or DFT) restricted = false } diff --git a/tests/li_scf_b3lyp_libxc/CMakeLists.txt b/tests/li_scf_b3lyp_libxc/CMakeLists.txt new file mode 100644 index 000000000..84592971d --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/CMakeLists.txt @@ -0,0 +1,11 @@ +if(ENABLE_MPI) + set(_li_scf_b3lyp_launcher "${MPIEXEC} ${MPIEXEC_NUMPROC_FLAG} 1") +endif() + +add_integration_test( + NAME "Li_SCF_B3LYP_Libxc" + LABELS "mrchem;libxc;Li_SCF_B3LYP_Libxc;dft;scf;energy" + COST 100 + LAUNCH_AGENT ${_li_scf_b3lyp_launcher} + INITIAL_GUESS ${CMAKE_CURRENT_LIST_DIR}/initial_guess + ) diff --git a/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.bas b/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.bas new file mode 100644 index 000000000..8297ac87f --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.bas @@ -0,0 +1,15 @@ +Gaussian basis 3-21G + 1 + 3. 1 2 1 1 +Li 0.0000000000 0.0000000000 0.0000000000 + 6 3 + 36.83820000 0.06966860 0.00000000 0.00000000 + 5.48172000 0.38134600 0.00000000 0.00000000 + 1.11327000 0.68170200 0.00000000 0.00000000 + 0.54020500 0.00000000 -0.26312700 0.00000000 + 0.10225500 0.00000000 1.14339000 0.00000000 + 0.02856500 0.00000000 0.00000000 1.00000000 + 3 2 + 0.54020500 0.16154600 0.00000000 + 0.10225500 0.91566300 0.00000000 + 0.02856500 0.00000000 1.00000000 diff --git a/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.moa b/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.moa new file mode 100644 index 000000000..58c1be2b3 --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.moa @@ -0,0 +1,82 @@ + 9 + -0.991216204923366 + -0.063168252106510 + 0.024611293917887 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.196235798467558 + 0.375192799624670 + 0.691835158906554 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.190262650474859 + 0.032894522967651 + 0.004820509091803 + -0.861003213282631 + 0.148858905853881 + 0.021814443388345 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.026206781579140 + -0.131161412253871 + -0.139337904569163 + -0.118594601163812 + -0.593550006411759 + -0.630551415460128 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.020457065978207 + 0.137912281574893 + -0.133667027187944 + 0.092575182242079 + 0.624099986470084 + -0.604888766300268 + -0.106844262187878 + 1.571170796830403 + -1.450216603570047 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.000000000000000 + -0.195731337542868 + 1.206022767761044 + -0.000000000000000 + 0.140221601585938 + -0.863992685931151 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 1.221802632528810 + -0.000000000000000 + 0.000000000000000 + -0.875297354556639 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 1.206022767761045 + 0.195731337542868 + 0.000000000000000 + -0.863992685931151 + -0.140221601585938 diff --git a/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.mob b/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.mob new file mode 100644 index 000000000..5aabab662 --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/initial_guess/mrchem.mob @@ -0,0 +1,82 @@ + 9 + -0.990741850222882 + -0.065977271109822 + 0.025495307561105 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.160323211242760 + -0.072975270859988 + 1.065158674860416 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.032062837182101 + 0.000000000000000 + 0.000000000000000 + -0.980791800359539 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.032062837182101 + 0.000000000000000 + 0.000000000000000 + -0.980791800359539 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.032062837182101 + 0.000000000000000 + 0.000000000000000 + -0.980791800359539 + 0.000000000000000 + -0.158621574885329 + 1.613585804138375 + -1.202978298565597 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -1.236559262809938 + 0.000000000000000 + 0.000000000000000 + 0.753760094669861 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -1.236559262809938 + 0.000000000000000 + 0.000000000000000 + 0.753760094669861 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -1.236559262809938 + 0.000000000000000 + 0.000000000000000 + 0.753760094669861 diff --git a/tests/li_scf_b3lyp_libxc/li.inp b/tests/li_scf_b3lyp_libxc/li.inp new file mode 100644 index 000000000..56dd12f76 --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/li.inp @@ -0,0 +1,32 @@ +# vim:syntax=sh: + +world_prec = 1.0e-3 # Overall relative precision +world_size = 5 # Size of simulation box 2^n + +MPI { + numerically_exact = true # Guarantee identical results in MPI +} + +Molecule { +multiplicity = 2 +$coords +Li 0.0 0.0 0.0 +$end +} + +WaveFunction { + method = DFT # Wave function method (HF or DFT) + restricted = false +} + +SCF { + kain = 3 # Length of KAIN iterative history + max_iter = 5 + orbital_thrs = 2.0e-2 + guess_type = GTO # Type of initial guess (none, gto, mw) +} + +DFT { + xc_library = libxc + functionals = b3lyp +} \ No newline at end of file diff --git a/tests/li_scf_b3lyp_libxc/li.out b/tests/li_scf_b3lyp_libxc/li.out new file mode 100644 index 000000000..fb7eb67ac --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/li.out @@ -0,0 +1,317 @@ + + +*************************************************************************** +*** *** +*** *** +*** __ __ ____ ____ _ *** +*** | \/ | _ \ / ___| |__ ___ _ __ ___ *** +*** | |\/| | |_) | | | '_ \ / _ \ '_ ` _ \ *** +*** | | | | _ <| |___| | | | __/ | | | | | *** +*** |_| |_|_| \_\\____|_| |_|\___|_| |_| |_| *** +*** *** +*** VERSION 1.2.0-alpha *** +*** *** +*** Git branch libxc_testbranch *** +*** Git commit hash 2bb9e2fa65c0876a30c8-dirty *** +*** Git commit author ylvao *** +*** Git commit date Thu Jan 29 12:29:18 2026 +0100 *** +*** *** +*** Contact: luca.frediani@uit.no *** +*** *** +*** Radovan Bast Magnar Bjorgve *** +*** Roberto Di Remigio Antoine Durdek *** +*** Luca Frediani Gabriel Gerez *** +*** Stig Rune Jensen Jonas Juselius *** +*** Rune Monstad Peter Wind *** +*** *** +*************************************************************************** + +--------------------------------------------------------------------------- + + MPI processes : (no bank) 1 + OpenMP threads : 1 + Total cores : 1 + +--------------------------------------------------------------------------- + + +--------------------------------------------------------------------------- + + MRCPP version : 1.6.0-alpha + Git branch : HEAD + Git commit hash : cec3cc96d578cdf77e8e + Git commit author : Niklas Göllmann + Git commit date : Tue Jan 6 16:29:30 2026 +0100 + + Linear algebra : EIGEN v3.4.0 + Parallelization : NONE + +--------------------------------------------------------------------------- + + + +=========================================================================== + MultiResolution Analysis +--------------------------------------------------------------------------- + polynomial order : 5 + polynomial type : Interpolating +--------------------------------------------------------------------------- + total boxes : 8 + boxes : [ 2 2 2 ] + unit lengths : [ 16.00000 16.00000 16.00000 ] + scaling factor : [ 1.00000 1.00000 1.00000 ] + lower bounds : [ -16.00000 -16.00000 -16.00000 ] + upper bounds : [ 16.00000 16.00000 16.00000 ] + total length : [ 32.00000 32.00000 32.00000 ] +=========================================================================== + + + +*************************************************************************** +*** *** +*** Initializing Molecule *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 2 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 Li : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +Name used in MapFunctionalName: B3LYP +Functional number: 0: b3lyp + + +*************************************************************************** +*** *** +*** XC Functional *** +*** *** +*************************************************************************** + +Using Libxc (version 7.0.0) to evaluate density functionals. Libxc is free +software. It is distributed under the Mozilla Public License, version 2.0. +For more information, please check the Libxc manual. You should cite + +S. Lehtola, C. Steigemann, M. J.T. Oliveira, and M. A.L. Marques., +SoftwareX 7, 1–5 (2018) DOI: 10.1016/j.softx.2017.11.002 + +when reporting the results of your calculation in a scientific article. + +Libxc functionals used in this calculation: + - hyb_gga_xc_b3lyp5 (ID 475): B3LYP with VWN functional 5 instead of RPA + * P. J. Stephens, F. J. Devlin, C. F. Chabalowski, and M. J. Frisch., J. Phys. Chem. 98, 11623 (1994) (DOI:10.1021/j100096a001) + +*************************************************************************** +*** *** +*** Computing Initial Guess Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial orbitals + Method : Project GTO molecular orbitals + Precision : 1.00000e-03 + Screening : 1.20000e+01 StdDev + Restricted : False + MO alpha file : initial_guess/mrchem.moa + MO beta file : initial_guess/mrchem.mob +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Orbitals +--------------------------------------------------------------------------- + Alpha electrons : 2 + Beta electrons : 1 + Total electrons : 3 +--------------------------------------------------------------------------- + n Occ Spin : Norm +--------------------------------------------------------------------------- + 0 1 a : 9.999999000474e-01 + 1 1 a : 1.000000190643e+00 + 2 1 b : 9.999998992098e-01 +=========================================================================== + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial energy + Method : DFT + XC Library : LibXC + Relativity : None + Environment : None + External fields : None + Precision : 1.00000e-03 + Density cutoff : 1.00000e-11 + Localization : Off +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +using libxc +=========================================================================== + Molecular Energy (initial) +--------------------------------------------------------------------------- + Kinetic energy : (au) 7.370997444552 + E-N energy : (au) -17.035029696825 + Coulomb energy : (au) 4.061916632919 + Exchange energy : (au) -0.355814809951 + X-C energy : (au) -1.469894461875 + N-N energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Electronic energy : (au) -7.427824891180 + Nuclear energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Total energy : (au) -7.427824891180e+00 + : (kcal/mol) -4.661030490897e+03 + : (kJ/mol) -1.950175157391e+04 + : (eV) -2.021214120815e+02 +=========================================================================== + + +=========================================================================== + Orbital Energies (initial) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 1 a : (au) -2.003629685543 + 1 1 a : (au) -0.131615364289 + 2 1 b : (au) -1.996727650980 +--------------------------------------------------------------------------- + Sum occupied : (au) -4.131972700812 +=========================================================================== + + + +*************************************************************************** +*** *** +*** Computing Ground State Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Optimize ground state orbitals + Method : DFT + XC library : LibXC + Relativity : None + Environment : None + External fields : None + Checkpointing : Off + Max iterations : 5 + KAIN solver : 3 + Localization : Off + Diagonalization : First two iterations + Start precision : 1.00000e-03 + Final precision : 1.00000e-03 + Helmholtz precision : Dynamic + Energy threshold : Off + Orbital threshold : 2.00000e-02 +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Iter MO residual Total energy Update +--------------------------------------------------------------------------- + 0 1.732051e+00 -7.427824891180 -7.427825e+00 + 1 6.957260e-02 -7.480527794641 -5.270290e-02 + 2 1.726775e-02 -7.481513818702 -9.860241e-04 +--------------------------------------------------------------------------- + SCF converged in 2 iterations! +=========================================================================== + + + +*************************************************************************** +*** *** +*** Printing Molecular Properties *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 2 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 Li : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +=========================================================================== + Molecular Energy (final) +--------------------------------------------------------------------------- + Kinetic energy : (au) 7.389078609209 + E-N energy : (au) -17.101945571040 + Coulomb energy : (au) 4.052896310158 + Exchange energy : (au) -0.354070574477 + X-C energy : (au) -1.467472592552 + N-N energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Electronic energy : (au) -7.481513818702 + Nuclear energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Total energy : (au) -7.481513818702e+00 + : (kcal/mol) -4.694720801569e+03 + : (kJ/mol) -1.964271183377e+04 + : (eV) -2.035823622254e+02 +=========================================================================== + + +=========================================================================== + Orbital Energies (final) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 1 a : (au) -2.033410667679 + 1 1 a : (au) -0.131567486907 + 2 1 b : (au) -2.024795910257 +--------------------------------------------------------------------------- + Sum occupied : (au) -4.189774064844 +=========================================================================== + + +=========================================================================== + Dipole Moment (dip-1) +--------------------------------------------------------------------------- + r_O : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Electronic vector : -0.000000 -0.000000 -0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Nuclear vector : 0.000000 0.000000 0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Total vector : 0.000000 0.000000 0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** Exiting MRChem *** +*** *** +*** Wall time : 0h 0m 25s *** +*** *** +*************************************************************************** + + diff --git a/tests/li_scf_b3lyp_libxc/reference/li.json b/tests/li_scf_b3lyp_libxc/reference/li.json new file mode 100644 index 000000000..9a5d922c3 --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/reference/li.json @@ -0,0 +1,336 @@ +{ + "input": { + "constants": { + "N_a": 6.02214076e+23, + "angstrom2bohrs": 1.8897261246257702, + "boltzmann_constant": 1.380649e-23, + "dipmom_au2debye": 2.5417464739297717, + "e0": 8.8541878128e-12, + "electron_g_factor": -2.00231930436256, + "elementary_charge": 1.602176634e-19, + "fine_structure_constant": 0.0072973525693, + "hartree2ev": 27.211386245988, + "hartree2kcalmol": 627.5094740630558, + "hartree2kjmol": 2625.4996394798254, + "hartree2simagnetizability": 78.9451185, + "hartree2wavenumbers": 219474.6313632, + "light_speed": 137.035999084, + "meter2bohr": 18897261246.2577 + }, + "geom_opt": { + "init_step_size": -0.5, + "max_force_component": 0.005, + "max_history_length": 10, + "max_iter": 100, + "minimal_step_size": 0.01, + "run": false, + "subspace_tolerance": 0.001, + "use_previous_guess": false + }, + "molecule": { + "charge": 0, + "coords": [ + { + "atom": "li", + "r_rms": 4.0992133976e-05, + "xyz": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "multiplicity": 2 + }, + "mpi": { + "bank_per_node": -1, + "bank_size": -1, + "numerically_exact": true, + "omp_threads": -1, + "shared_memory_size": 10000, + "use_omp_num_threads": false + }, + "mra": { + "basis_order": 5, + "basis_type": "interpolating", + "boxes": [ + 2, + 2, + 2 + ], + "corner": [ + -1, + -1, + -1 + ], + "max_scale": 20, + "min_scale": -4 + }, + "printer": { + "file_name": "li", + "print_constants": false, + "print_level": 0, + "print_mpi": false, + "print_prec": 6, + "print_width": 75 + }, + "rsp_calculations": {}, + "scf_calculation": { + "fock_operator": { + "coulomb_operator": { + "poisson_prec": 0.001, + "shared_memory": false + }, + "exchange_operator": { + "exchange_prec": -1.0, + "poisson_prec": 0.001 + }, + "kinetic_operator": { + "derivative": "abgv_55" + }, + "nuclear_operator": { + "nuclear_model": "point_like", + "proj_prec": 0.001, + "shared_memory": false, + "smooth_prec": 0.001 + }, + "xc_library": [ + "xcfun" + ], + "xc_operator": { + "shared_memory": false, + "xc_functional": { + "cutoff": 1e-11, + "functionals": [ + { + "coef": 1.0, + "name": "b3lyp" + } + ], + "spin": true + } + } + }, + "initial_guess": { + "cutoff": 1e-11, + "environment": "None", + "external_field": "None", + "file_CUBE_a": "cube_vectors/CUBE_a_vector.json", + "file_CUBE_b": "cube_vectors/CUBE_b_vector.json", + "file_CUBE_p": "cube_vectors/CUBE_p_vector.json", + "file_basis": "initial_guess/mrchem.bas", + "file_chk": "checkpoint/phi_scf", + "file_gto_a": "initial_guess/mrchem.moa", + "file_gto_b": "initial_guess/mrchem.mob", + "file_gto_p": "initial_guess/mrchem.mop", + "file_phi_a": "initial_guess/phi_a_scf", + "file_phi_b": "initial_guess/phi_b_scf", + "file_phi_p": "initial_guess/phi_p_scf", + "localize": false, + "method": "DFT", + "prec": 0.001, + "relativity": "None", + "restricted": false, + "rotate": true, + "screen": 12.0, + "type": "gto", + "xc_library": "xcfun", + "zeta": 0 + }, + "properties": { + "dipole_moment": { + "dip-1": { + "operator": "h_e_dip", + "precision": 0.001, + "r_O": [ + 0.0, + 0.0, + 0.0 + ] + } + } + }, + "scf_solver": { + "checkpoint": false, + "energy_thrs": -1.0, + "environment": "None", + "external_field": "None", + "file_chk": "checkpoint/phi_scf", + "final_prec": 0.001, + "helmholtz_prec": -1.0, + "kain": 3, + "localize": false, + "max_iter": 5, + "method": "DFT", + "orbital_thrs": 0.02, + "relativity": "None", + "rotation": 0, + "start_prec": 0.001, + "xc_library": "xcfun" + } + }, + "schema_name": "mrchem_input", + "schema_version": 1 + }, + "output": { + "properties": { + "center_of_mass": [ + 0.0, + 0.0, + 0.0 + ], + "charge": 0, + "dipole_moment": { + "dip-1": { + "magnitude": 6.63546194208962e-14, + "r_O": [ + 0.0, + 0.0, + 0.0 + ], + "vector": [ + 0.0, + 0.0, + 0.0 + ], + "vector_el": [ + 0.0, + 0.0, + 0.0 + ], + "vector_nuc": [ + 0.0, + 0.0, + 0.0 + ] + } + }, + "geometry": [ + { + "symbol": "Li", + "xyz": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "multiplicity": 2, + "orbital_energies": { + "energy": [ + -2.0334106780827024, + -0.13156748738892196, + -2.0247959175562977 + ], + "occupation": [ + 1.0, + 1.0, + 1.0 + ], + "spin": [ + "a", + "a", + "b" + ], + "sum_occupied": -4.189774083027922 + }, + "scf_energy": { + "E_ee": 4.052896289844306, + "E_eext": 0.0, + "E_el": -7.481513824145804, + "E_en": -17.1019455092072, + "E_kin": 7.389078559562881, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.481513824145804, + "E_x": -0.3540705723345492, + "E_xc": -1.4674725920112404, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + } + }, + "provenance": { + "creator": "MRChem", + "mpi_processes": 1, + "nthreads": 1, + "routine": "mrchem.x", + "total_cores": 1, + "version": "1.2.0-alpha" + }, + "rsp_calculations": null, + "scf_calculation": { + "initial_energy": { + "E_ee": 4.061916632919257, + "E_eext": 0.0, + "E_el": -7.4278248937059335, + "E_en": -17.035029696824992, + "E_kin": 7.370997444552398, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.4278248937059335, + "E_x": -0.35581480995128667, + "E_xc": -1.4698944644013088, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "scf_solver": { + "converged": true, + "cycles": [ + { + "energy_terms": { + "E_ee": 4.0372088079675965, + "E_eext": 0.0, + "E_el": -7.480527798944603, + "E_en": -17.01787840416016, + "E_kin": 7.314949642048837, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.480527798944603, + "E_x": -0.3529043821092998, + "E_xc": -1.4619034626915766, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -7.480527798944603, + "energy_update": 0.05270290523866983, + "mo_residual": 0.06957255189837697, + "wall_time": 21.121093435 + }, + { + "energy_terms": { + "E_ee": 4.052896289844306, + "E_eext": 0.0, + "E_el": -7.481513824145804, + "E_en": -17.1019455092072, + "E_kin": 7.389078559562881, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.481513824145804, + "E_x": -0.3540705723345492, + "E_xc": -1.4674725920112404, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -7.481513824145804, + "energy_update": 0.0009860252012003912, + "mo_residual": 0.017267732458522325, + "wall_time": 13.443850889 + } + ], + "wall_time": 34.56734884 + }, + "success": true + }, + "schema_name": "mrchem_output", + "schema_version": 1, + "success": true + } +} diff --git a/tests/li_scf_b3lyp_libxc/reference/li.out b/tests/li_scf_b3lyp_libxc/reference/li.out new file mode 100644 index 000000000..49cbd0f65 --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/reference/li.out @@ -0,0 +1,315 @@ + + +*************************************************************************** +*** *** +*** *** +*** __ __ ____ ____ _ *** +*** | \/ | _ \ / ___| |__ ___ _ __ ___ *** +*** | |\/| | |_) | | | '_ \ / _ \ '_ ` _ \ *** +*** | | | | _ <| |___| | | | __/ | | | | | *** +*** |_| |_|_| \_\\____|_| |_|\___|_| |_| |_| *** +*** *** +*** VERSION 1.2.0-alpha *** +*** *** +*** Git branch libxc_testbranch *** +*** Git commit hash e6202c5abcba6d1f62f6-dirty *** +*** Git commit author ylvao *** +*** Git commit date Thu Jan 29 08:19:49 2026 +0100 *** +*** *** +*** Contact: luca.frediani@uit.no *** +*** *** +*** Radovan Bast Magnar Bjorgve *** +*** Roberto Di Remigio Antoine Durdek *** +*** Luca Frediani Gabriel Gerez *** +*** Stig Rune Jensen Jonas Juselius *** +*** Rune Monstad Peter Wind *** +*** *** +*************************************************************************** + +--------------------------------------------------------------------------- + + MPI processes : (no bank) 1 + OpenMP threads : 1 + Total cores : 1 + +--------------------------------------------------------------------------- + + +--------------------------------------------------------------------------- + + MRCPP version : 1.6.0-alpha + Git branch : HEAD + Git commit hash : cec3cc96d578cdf77e8e + Git commit author : Niklas Göllmann + Git commit date : Tue Jan 6 16:29:30 2026 +0100 + + Linear algebra : EIGEN v3.4.0 + Parallelization : NONE + +--------------------------------------------------------------------------- + + + +=========================================================================== + MultiResolution Analysis +--------------------------------------------------------------------------- + polynomial order : 5 + polynomial type : Interpolating +--------------------------------------------------------------------------- + total boxes : 8 + boxes : [ 2 2 2 ] + unit lengths : [ 16.00000 16.00000 16.00000 ] + scaling factor : [ 1.00000 1.00000 1.00000 ] + lower bounds : [ -16.00000 -16.00000 -16.00000 ] + upper bounds : [ 16.00000 16.00000 16.00000 ] + total length : [ 32.00000 32.00000 32.00000 ] +=========================================================================== + + + +*************************************************************************** +*** *** +*** Initializing Molecule *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 2 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 Li : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** XC Functional *** +*** *** +*************************************************************************** + +XCFun DFT library Copyright 2009-2020 Ulf Ekstrom and contributors. +See http://dftlibs.org/xcfun/ for more information. + +This is free software; see the source code for copying conditions. +There is ABSOLUTELY NO WARRANTY; not even for MERCHANTABILITY or +FITNESS FOR A PARTICULAR PURPOSE. For details see the documentation. +Scientific users of this library should cite +U. Ekstrom, L. Visscher, R. Bast, A. J. Thorvaldsen and K. Ruud; +J.Chem.Theor.Comp. 2010, DOI: 10.1021/ct100117s + +--------------------------------------------------------------------------- + + +*************************************************************************** +*** *** +*** Computing Initial Guess Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial orbitals + Method : Project GTO molecular orbitals + Precision : 1.00000e-03 + Screening : 1.20000e+01 StdDev + Restricted : False + MO alpha file : initial_guess/mrchem.moa + MO beta file : initial_guess/mrchem.mob +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Orbitals +--------------------------------------------------------------------------- + Alpha electrons : 2 + Beta electrons : 1 + Total electrons : 3 +--------------------------------------------------------------------------- + n Occ Spin : Norm +--------------------------------------------------------------------------- + 0 1 a : 9.999999000474e-01 + 1 1 a : 1.000000190643e+00 + 2 1 b : 9.999998992098e-01 +=========================================================================== + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial energy + Method : DFT + XC Library : XCFun + Relativity : None + Environment : None + External fields : None + Precision : 1.00000e-03 + Density cutoff : 1.00000e-11 + Localization : Off +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +using xcfun +=========================================================================== + Molecular Energy (initial) +--------------------------------------------------------------------------- + Kinetic energy : (au) 7.370997444552 + E-N energy : (au) -17.035029696825 + Coulomb energy : (au) 4.061916632919 + Exchange energy : (au) -0.355814809951 + X-C energy : (au) -1.469894464401 + N-N energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Electronic energy : (au) -7.427824893706 + Nuclear energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Total energy : (au) -7.427824893706e+00 + : (kcal/mol) -4.661030492482e+03 + : (kJ/mol) -1.950175158054e+04 + : (eV) -2.021214121502e+02 +=========================================================================== + + +=========================================================================== + Orbital Energies (initial) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 1 a : (au) -2.003629685543 + 1 1 a : (au) -0.131615364274 + 2 1 b : (au) -1.996727661485 +--------------------------------------------------------------------------- + Sum occupied : (au) -4.131972711302 +=========================================================================== + + + +*************************************************************************** +*** *** +*** Computing Ground State Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Optimize ground state orbitals + Method : DFT + XC library : XCFun + Relativity : None + Environment : None + External fields : None + Checkpointing : Off + Max iterations : 5 + KAIN solver : 3 + Localization : Off + Diagonalization : First two iterations + Start precision : 1.00000e-03 + Final precision : 1.00000e-03 + Helmholtz precision : Dynamic + Energy threshold : Off + Orbital threshold : 2.00000e-02 +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Iter MO residual Total energy Update +--------------------------------------------------------------------------- + 0 1.732051e+00 -7.427824893706 -7.427825e+00 + 1 6.957255e-02 -7.480527798945 -5.270291e-02 + 2 1.726773e-02 -7.481513824146 -9.860252e-04 +--------------------------------------------------------------------------- + SCF converged in 2 iterations! +=========================================================================== + + + +*************************************************************************** +*** *** +*** Printing Molecular Properties *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 2 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 Li : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +=========================================================================== + Molecular Energy (final) +--------------------------------------------------------------------------- + Kinetic energy : (au) 7.389078559563 + E-N energy : (au) -17.101945509207 + Coulomb energy : (au) 4.052896289844 + Exchange energy : (au) -0.354070572335 + X-C energy : (au) -1.467472592011 + N-N energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Electronic energy : (au) -7.481513824146 + Nuclear energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Total energy : (au) -7.481513824146e+00 + : (kcal/mol) -4.694720804985e+03 + : (kJ/mol) -1.964271184806e+04 + : (eV) -2.035823623735e+02 +=========================================================================== + + +=========================================================================== + Orbital Energies (final) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 1 a : (au) -2.033410678083 + 1 1 a : (au) -0.131567487389 + 2 1 b : (au) -2.024795917556 +--------------------------------------------------------------------------- + Sum occupied : (au) -4.189774083028 +=========================================================================== + + +=========================================================================== + Dipole Moment (dip-1) +--------------------------------------------------------------------------- + r_O : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Electronic vector : -0.000000 -0.000000 -0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Nuclear vector : 0.000000 0.000000 0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Total vector : 0.000000 0.000000 0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** Exiting MRChem *** +*** *** +*** Wall time : 0h 0m 43s *** +*** *** +*************************************************************************** + + diff --git a/tests/li_scf_b3lyp_libxc/test b/tests/li_scf_b3lyp_libxc/test new file mode 100755 index 000000000..5aff74384 --- /dev/null +++ b/tests/li_scf_b3lyp_libxc/test @@ -0,0 +1,26 @@ +#!/usr/bin/env python3 + +import sys +from pathlib import Path + +sys.path.append(str(Path(__file__).resolve().parents[1])) + +from tester import * # isort:skip + +options = script_cli() + +filters = { + SUM_OCCUPIED: rel_tolerance(1.0e-7), + E_KIN: rel_tolerance(1.0e-7), + E_EN: rel_tolerance(1.0e-7), + E_EE: rel_tolerance(1.0e-7), + E_X: rel_tolerance(1.0e-7), + E_XC: rel_tolerance(1.0e-7), + E_EEXT: rel_tolerance(1.0e-7), + E_NEXT: rel_tolerance(1.0e-7), + E_EL: rel_tolerance(1.0e-7), +} + +ierr = run(options, input_file="li", filters=filters) + +sys.exit(ierr) diff --git a/tests/li_scf_lda_libxc/li.inp b/tests/li_scf_lda_libxc/li.inp new file mode 100644 index 000000000..0acc3b28a --- /dev/null +++ b/tests/li_scf_lda_libxc/li.inp @@ -0,0 +1,31 @@ +# vim:syntax=sh: + +world_prec = 1.0e-3 # Overall relative precision +world_size = 5 # Size of simulation box 2^n + +MPI { + numerically_exact = true # Guarantee identical results in MPI +} + +Molecule { +multiplicity = 2 +$coords +Li 0.0 0.0 0.0 +$end +} + +WaveFunction { + method = SVWN5 # Wave function method (HF or DFT) + restricted = false +} + +SCF { + kain = 3 # Length of KAIN iterative history + max_iter = 5 + orbital_thrs = 2.0e-2 + guess_type = GTO # Type of initial guess (none, gto, mw) +} + +DFT { + xc_library = libxc +} \ No newline at end of file diff --git a/tests/li_scf_svwn5_libxc/CMakeLists.txt b/tests/li_scf_svwn5_libxc/CMakeLists.txt new file mode 100644 index 000000000..8beb4af56 --- /dev/null +++ b/tests/li_scf_svwn5_libxc/CMakeLists.txt @@ -0,0 +1,11 @@ +if(ENABLE_MPI) + set(_li_scf_svwn5_libxc_launcher "${MPIEXEC} ${MPIEXEC_NUMPROC_FLAG} 1") +endif() + +add_integration_test( + NAME "Li_SCF_SVWN5_Libxc" + LABELS "mrchem;li_scf_svwn5_libxc;Li_SCF_SVWN5_Libxc;libxc;dft;scf;energy;" + COST 100 + LAUNCH_AGENT ${_li_scf_svwn5_libxc_launcher} + INITIAL_GUESS ${CMAKE_CURRENT_LIST_DIR}/initial_guess + ) diff --git a/tests/li_scf_svwn5_libxc/initial_guess/mrchem.bas b/tests/li_scf_svwn5_libxc/initial_guess/mrchem.bas new file mode 100644 index 000000000..8297ac87f --- /dev/null +++ b/tests/li_scf_svwn5_libxc/initial_guess/mrchem.bas @@ -0,0 +1,15 @@ +Gaussian basis 3-21G + 1 + 3. 1 2 1 1 +Li 0.0000000000 0.0000000000 0.0000000000 + 6 3 + 36.83820000 0.06966860 0.00000000 0.00000000 + 5.48172000 0.38134600 0.00000000 0.00000000 + 1.11327000 0.68170200 0.00000000 0.00000000 + 0.54020500 0.00000000 -0.26312700 0.00000000 + 0.10225500 0.00000000 1.14339000 0.00000000 + 0.02856500 0.00000000 0.00000000 1.00000000 + 3 2 + 0.54020500 0.16154600 0.00000000 + 0.10225500 0.91566300 0.00000000 + 0.02856500 0.00000000 1.00000000 diff --git a/tests/li_scf_svwn5_libxc/initial_guess/mrchem.moa b/tests/li_scf_svwn5_libxc/initial_guess/mrchem.moa new file mode 100644 index 000000000..58c1be2b3 --- /dev/null +++ b/tests/li_scf_svwn5_libxc/initial_guess/mrchem.moa @@ -0,0 +1,82 @@ + 9 + -0.991216204923366 + -0.063168252106510 + 0.024611293917887 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.196235798467558 + 0.375192799624670 + 0.691835158906554 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.190262650474859 + 0.032894522967651 + 0.004820509091803 + -0.861003213282631 + 0.148858905853881 + 0.021814443388345 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.026206781579140 + -0.131161412253871 + -0.139337904569163 + -0.118594601163812 + -0.593550006411759 + -0.630551415460128 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.020457065978207 + 0.137912281574893 + -0.133667027187944 + 0.092575182242079 + 0.624099986470084 + -0.604888766300268 + -0.106844262187878 + 1.571170796830403 + -1.450216603570047 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.000000000000000 + -0.195731337542868 + 1.206022767761044 + -0.000000000000000 + 0.140221601585938 + -0.863992685931151 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 1.221802632528810 + -0.000000000000000 + 0.000000000000000 + -0.875297354556639 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 1.206022767761045 + 0.195731337542868 + 0.000000000000000 + -0.863992685931151 + -0.140221601585938 diff --git a/tests/li_scf_svwn5_libxc/initial_guess/mrchem.mob b/tests/li_scf_svwn5_libxc/initial_guess/mrchem.mob new file mode 100644 index 000000000..5aabab662 --- /dev/null +++ b/tests/li_scf_svwn5_libxc/initial_guess/mrchem.mob @@ -0,0 +1,82 @@ + 9 + -0.990741850222882 + -0.065977271109822 + 0.025495307561105 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.160323211242760 + -0.072975270859988 + 1.065158674860416 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.032062837182101 + 0.000000000000000 + 0.000000000000000 + -0.980791800359539 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.032062837182101 + 0.000000000000000 + 0.000000000000000 + -0.980791800359539 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -0.032062837182101 + 0.000000000000000 + 0.000000000000000 + -0.980791800359539 + 0.000000000000000 + -0.158621574885329 + 1.613585804138375 + -1.202978298565597 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -1.236559262809938 + 0.000000000000000 + 0.000000000000000 + 0.753760094669861 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -1.236559262809938 + 0.000000000000000 + 0.000000000000000 + 0.753760094669861 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + 0.000000000000000 + -1.236559262809938 + 0.000000000000000 + 0.000000000000000 + 0.753760094669861 diff --git a/tests/li_scf_svwn5_libxc/li.inp b/tests/li_scf_svwn5_libxc/li.inp new file mode 100644 index 000000000..65f2ec247 --- /dev/null +++ b/tests/li_scf_svwn5_libxc/li.inp @@ -0,0 +1,31 @@ +# vim:syntax=sh: + +world_prec = 1.0e-3 # Overall relative precision +world_size = 5 # Size of simulation box 2^n + +MPI { + numerically_exact = true # Guarantee identical results in MPI +} + +Molecule { +multiplicity = 2 +$coords +Li 0.0 0.0 0.0 +$end +} + +WaveFunction { + method = svwn5 # Wave function method (HF or DFT) + restricted = false +} + +SCF { + kain = 3 # Length of KAIN iterative history + max_iter = 5 + orbital_thrs = 2.0e-2 + guess_type = GTO # Type of initial guess (none, gto, mw) +} + +DFT { + xc_library = libxc +} \ No newline at end of file diff --git a/tests/li_scf_svwn5_libxc/reference/li.json b/tests/li_scf_svwn5_libxc/reference/li.json new file mode 100644 index 000000000..b4ee88b84 --- /dev/null +++ b/tests/li_scf_svwn5_libxc/reference/li.json @@ -0,0 +1,336 @@ +{ + "input": { + "constants": { + "N_a": 6.02214076e+23, + "angstrom2bohrs": 1.8897261246257702, + "boltzmann_constant": 1.380649e-23, + "dipmom_au2debye": 2.5417464739297717, + "e0": 8.8541878128e-12, + "electron_g_factor": -2.00231930436256, + "elementary_charge": 1.602176634e-19, + "fine_structure_constant": 0.0072973525693, + "hartree2ev": 27.211386245988, + "hartree2kcalmol": 627.5094740630558, + "hartree2kjmol": 2625.4996394798254, + "hartree2simagnetizability": 78.9451185, + "hartree2wavenumbers": 219474.6313632, + "light_speed": 137.035999084, + "meter2bohr": 18897261246.2577 + }, + "geom_opt": { + "init_step_size": -0.5, + "max_force_component": 0.005, + "max_history_length": 10, + "max_iter": 100, + "minimal_step_size": 0.01, + "run": false, + "subspace_tolerance": 0.001, + "use_previous_guess": false + }, + "molecule": { + "charge": 0, + "coords": [ + { + "atom": "li", + "r_rms": 4.0992133976e-05, + "xyz": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "multiplicity": 2 + }, + "mpi": { + "bank_per_node": -1, + "bank_size": -1, + "numerically_exact": true, + "omp_threads": -1, + "shared_memory_size": 10000, + "use_omp_num_threads": false + }, + "mra": { + "basis_order": 5, + "basis_type": "interpolating", + "boxes": [ + 2, + 2, + 2 + ], + "corner": [ + -1, + -1, + -1 + ], + "max_scale": 20, + "min_scale": -4 + }, + "printer": { + "file_name": "li", + "print_constants": false, + "print_level": 0, + "print_mpi": false, + "print_prec": 6, + "print_width": 75 + }, + "rsp_calculations": {}, + "scf_calculation": { + "fock_operator": { + "coulomb_operator": { + "poisson_prec": 0.001, + "shared_memory": false + }, + "exchange_operator": { + "exchange_prec": -1.0, + "poisson_prec": 0.001 + }, + "kinetic_operator": { + "derivative": "abgv_55" + }, + "nuclear_operator": { + "nuclear_model": "point_like", + "proj_prec": 0.001, + "shared_memory": false, + "smooth_prec": 0.001 + }, + "xc_library": [ + "xcfun" + ], + "xc_operator": { + "shared_memory": false, + "xc_functional": { + "cutoff": 1e-11, + "functionals": [ + { + "coef": 1.0, + "name": "svwn5" + } + ], + "spin": true + } + } + }, + "initial_guess": { + "cutoff": 1e-11, + "environment": "None", + "external_field": "None", + "file_CUBE_a": "cube_vectors/CUBE_a_vector.json", + "file_CUBE_b": "cube_vectors/CUBE_b_vector.json", + "file_CUBE_p": "cube_vectors/CUBE_p_vector.json", + "file_basis": "initial_guess/mrchem.bas", + "file_chk": "checkpoint/phi_scf", + "file_gto_a": "initial_guess/mrchem.moa", + "file_gto_b": "initial_guess/mrchem.mob", + "file_gto_p": "initial_guess/mrchem.mop", + "file_phi_a": "initial_guess/phi_a_scf", + "file_phi_b": "initial_guess/phi_b_scf", + "file_phi_p": "initial_guess/phi_p_scf", + "localize": false, + "method": "DFT (SVWN5)", + "prec": 0.001, + "relativity": "None", + "restricted": false, + "rotate": true, + "screen": 12.0, + "type": "gto", + "xc_library": "xcfun", + "zeta": 0 + }, + "properties": { + "dipole_moment": { + "dip-1": { + "operator": "h_e_dip", + "precision": 0.001, + "r_O": [ + 0.0, + 0.0, + 0.0 + ] + } + } + }, + "scf_solver": { + "checkpoint": false, + "energy_thrs": -1.0, + "environment": "None", + "external_field": "None", + "file_chk": "checkpoint/phi_scf", + "final_prec": 0.001, + "helmholtz_prec": -1.0, + "kain": 3, + "localize": false, + "max_iter": 5, + "method": "DFT (SVWN5)", + "orbital_thrs": 0.02, + "relativity": "None", + "rotation": 0, + "start_prec": 0.001, + "xc_library": "xcfun" + } + }, + "schema_name": "mrchem_input", + "schema_version": 1 + }, + "output": { + "properties": { + "center_of_mass": [ + 0.0, + 0.0, + 0.0 + ], + "charge": 0, + "dipole_moment": { + "dip-1": { + "magnitude": 7.347868054993544e-14, + "r_O": [ + 0.0, + 0.0, + 0.0 + ], + "vector": [ + 0.0, + 0.0, + 0.0 + ], + "vector_el": [ + 0.0, + 0.0, + 0.0 + ], + "vector_nuc": [ + 0.0, + 0.0, + 0.0 + ] + } + }, + "geometry": [ + { + "symbol": "Li", + "xyz": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "multiplicity": 2, + "orbital_energies": { + "energy": [ + -1.8772810067255639, + -0.1166720344906071, + -1.8696802475837315 + ], + "occupation": [ + 1.0, + 1.0, + 1.0 + ], + "spin": [ + "a", + "a", + "b" + ], + "sum_occupied": -3.8636332887999023 + }, + "scf_energy": { + "E_ee": 4.0022303119866836, + "E_eext": 0.0, + "E_el": -7.342983769075873, + "E_en": -16.902251628428512, + "E_kin": 7.219654939850632, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.342983769075873, + "E_x": 0.0, + "E_xc": -1.6626173924846748, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + } + }, + "provenance": { + "creator": "MRChem", + "mpi_processes": 1, + "nthreads": 64, + "routine": "mrchem.x", + "total_cores": 64, + "version": "1.2.0-alpha" + }, + "rsp_calculations": null, + "scf_calculation": { + "initial_energy": { + "E_ee": 4.061916632919257, + "E_eext": 0.0, + "E_el": -7.288838366461443, + "E_en": -17.035029696824992, + "E_kin": 7.370997444552398, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.288838366461443, + "E_x": 0.0, + "E_xc": -1.6867227471081052, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "scf_solver": { + "converged": true, + "cycles": [ + { + "energy_terms": { + "E_ee": 3.9896186379178404, + "E_eext": 0.0, + "E_el": -7.34236754400743, + "E_en": -16.837854221195187, + "E_kin": 7.163926304298668, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.34236754400743, + "E_x": 0.0, + "E_xc": -1.6580582650287512, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -7.34236754400743, + "energy_update": 0.05352917754598696, + "mo_residual": 0.07476458299670477, + "wall_time": 5.493513328 + }, + { + "energy_terms": { + "E_ee": 4.0022303119866836, + "E_eext": 0.0, + "E_el": -7.342983769075873, + "E_en": -16.902251628428512, + "E_kin": 7.219654939850632, + "E_next": 0.0, + "E_nn": 0.0, + "E_nuc": 0.0, + "E_tot": -7.342983769075873, + "E_x": 0.0, + "E_xc": -1.6626173924846748, + "Er_el": 0.0, + "Er_nuc": 0.0, + "Er_tot": 0.0 + }, + "energy_total": -7.342983769075873, + "energy_update": 0.0006162250684429438, + "mo_residual": 0.016750742022743247, + "wall_time": 7.530144402 + } + ], + "wall_time": 13.02728183 + }, + "success": true + }, + "schema_name": "mrchem_output", + "schema_version": 1, + "success": true + } +} diff --git a/tests/li_scf_svwn5_libxc/reference/li.out b/tests/li_scf_svwn5_libxc/reference/li.out new file mode 100644 index 000000000..0962bdf48 --- /dev/null +++ b/tests/li_scf_svwn5_libxc/reference/li.out @@ -0,0 +1,314 @@ + + +*************************************************************************** +*** *** +*** *** +*** __ __ ____ ____ _ *** +*** | \/ | _ \ / ___| |__ ___ _ __ ___ *** +*** | |\/| | |_) | | | '_ \ / _ \ '_ ` _ \ *** +*** | | | | _ <| |___| | | | __/ | | | | | *** +*** |_| |_|_| \_\\____|_| |_|\___|_| |_| |_| *** +*** *** +*** VERSION 1.2.0-alpha *** +*** *** +*** Git branch libxc_testbranch *** +*** Git commit hash cdd0549cf55900c1b4d7 *** +*** Git commit author ylvao *** +*** Git commit date Mon Feb 2 15:24:40 2026 +0100 *** +*** *** +*** Contact: luca.frediani@uit.no *** +*** *** +*** Radovan Bast Magnar Bjorgve *** +*** Roberto Di Remigio Antoine Durdek *** +*** Luca Frediani Gabriel Gerez *** +*** Stig Rune Jensen Jonas Juselius *** +*** Rune Monstad Peter Wind *** +*** *** +*************************************************************************** + +--------------------------------------------------------------------------- + + MPI processes : (no bank) 1 + OpenMP threads : 64 + Total cores : 64 + +--------------------------------------------------------------------------- + + +--------------------------------------------------------------------------- + + MRCPP version : 1.6.0-alpha + Git branch : HEAD + Git commit hash : cec3cc96d578cdf77e8e + Git commit author : Niklas Göllmann + Git commit date : Tue Jan 6 16:29:30 2026 +0100 + + Linear algebra : EIGEN v3.4.0 + Parallelization : OpenMP (64 threads) + +--------------------------------------------------------------------------- + + + +=========================================================================== + MultiResolution Analysis +--------------------------------------------------------------------------- + polynomial order : 5 + polynomial type : Interpolating +--------------------------------------------------------------------------- + total boxes : 8 + boxes : [ 2 2 2 ] + unit lengths : [ 16.00000 16.00000 16.00000 ] + scaling factor : [ 1.00000 1.00000 1.00000 ] + lower bounds : [ -16.00000 -16.00000 -16.00000 ] + upper bounds : [ 16.00000 16.00000 16.00000 ] + total length : [ 32.00000 32.00000 32.00000 ] +=========================================================================== + + + +*************************************************************************** +*** *** +*** Initializing Molecule *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 2 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 Li : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** XC Functional *** +*** *** +*************************************************************************** + +XCFun DFT library Copyright 2009-2020 Ulf Ekstrom and contributors. +See http://dftlibs.org/xcfun/ for more information. + +This is free software; see the source code for copying conditions. +There is ABSOLUTELY NO WARRANTY; not even for MERCHANTABILITY or +FITNESS FOR A PARTICULAR PURPOSE. For details see the documentation. +Scientific users of this library should cite +U. Ekstrom, L. Visscher, R. Bast, A. J. Thorvaldsen and K. Ruud; +J.Chem.Theor.Comp. 2010, DOI: 10.1021/ct100117s + +--------------------------------------------------------------------------- + + +*************************************************************************** +*** *** +*** Computing Initial Guess Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial orbitals + Method : Project GTO molecular orbitals + Precision : 1.00000e-03 + Screening : 1.20000e+01 StdDev + Restricted : False + MO alpha file : initial_guess/mrchem.moa + MO beta file : initial_guess/mrchem.mob +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Orbitals +--------------------------------------------------------------------------- + Alpha electrons : 2 + Beta electrons : 1 + Total electrons : 3 +--------------------------------------------------------------------------- + n Occ Spin : Norm +--------------------------------------------------------------------------- + 0 1 a : 9.999999000474e-01 + 1 1 a : 1.000000190643e+00 + 2 1 b : 9.999998992098e-01 +=========================================================================== + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Compute initial energy + Method : DFT (SVWN5) + XC Library : XCFun + Relativity : None + Environment : None + External fields : None + Precision : 1.00000e-03 + Density cutoff : 1.00000e-11 + Localization : Off +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Molecular Energy (initial) +--------------------------------------------------------------------------- + Kinetic energy : (au) 7.370997444552 + E-N energy : (au) -17.035029696825 + Coulomb energy : (au) 4.061916632919 + Exchange energy : (au) 0.000000000000 + X-C energy : (au) -1.686722747108 + N-N energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Electronic energy : (au) -7.288838366461 + Nuclear energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Total energy : (au) -7.288838366461e+00 + : (kcal/mol) -4.573815129869e+03 + : (kJ/mol) -1.913684250337e+04 + : (eV) -1.983393960744e+02 +=========================================================================== + + +=========================================================================== + Orbital Energies (initial) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 1 a : (au) -1.824986094532 + 1 1 a : (au) -0.114890316102 + 2 1 b : (au) -1.817727169098 +--------------------------------------------------------------------------- + Sum occupied : (au) -3.757603579732 +=========================================================================== + + + +*************************************************************************** +*** *** +*** Computing Ground State Wavefunction *** +*** *** +*************************************************************************** + + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + Calculation : Optimize ground state orbitals + Method : DFT (SVWN5) + XC library : XCFun + Relativity : None + Environment : None + External fields : None + Checkpointing : Off + Max iterations : 5 + KAIN solver : 3 + Localization : Off + Diagonalization : First two iterations + Start precision : 1.00000e-03 + Final precision : 1.00000e-03 + Helmholtz precision : Dynamic + Energy threshold : Off + Orbital threshold : 2.00000e-02 +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + +=========================================================================== + Iter MO residual Total energy Update +--------------------------------------------------------------------------- + 0 1.732051e+00 -7.288838366461 -7.288838e+00 + 1 7.476458e-02 -7.342367544007 -5.352918e-02 + 2 1.675074e-02 -7.342983769076 -6.162251e-04 +--------------------------------------------------------------------------- + SCF converged in 2 iterations! +=========================================================================== + + + +*************************************************************************** +*** *** +*** Printing Molecular Properties *** +*** *** +*************************************************************************** + + +=========================================================================== + Molecule +--------------------------------------------------------------------------- + Charge : 0 + Multiplicity : 2 +--------------------------------------------------------------------------- + N Atom : x y z +--------------------------------------------------------------------------- + 0 Li : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Center of mass : 0.000000 0.000000 0.000000 +=========================================================================== + + +=========================================================================== + Molecular Energy (final) +--------------------------------------------------------------------------- + Kinetic energy : (au) 7.219654939851 + E-N energy : (au) -16.902251628429 + Coulomb energy : (au) 4.002230311987 + Exchange energy : (au) 0.000000000000 + X-C energy : (au) -1.662617392485 + N-N energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Electronic energy : (au) -7.342983769076 + Nuclear energy : (au) 0.000000000000 +--------------------------------------------------------------------------- + Total energy : (au) -7.342983769076e+00 + : (kcal/mol) -4.607791882986e+03 + : (kJ/mol) -1.927900123841e+04 + : (eV) -1.998127675383e+02 +=========================================================================== + + +=========================================================================== + Orbital Energies (final) +--------------------------------------------------------------------------- + n Occ Spin : Epsilon +--------------------------------------------------------------------------- + 0 1 a : (au) -1.877281006726 + 1 1 a : (au) -0.116672034491 + 2 1 b : (au) -1.869680247584 +--------------------------------------------------------------------------- + Sum occupied : (au) -3.863633288800 +=========================================================================== + + +=========================================================================== + Dipole Moment (dip-1) +--------------------------------------------------------------------------- + r_O : 0.000000 0.000000 0.000000 +--------------------------------------------------------------------------- + Electronic vector : -0.000000 -0.000000 -0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Nuclear vector : 0.000000 0.000000 0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +--------------------------------------------------------------------------- + Total vector : 0.000000 0.000000 0.000000 + Magnitude : (au) 0.000000 + : (Debye) 0.000000 +=========================================================================== + + + + +*************************************************************************** +*** *** +*** Exiting MRChem *** +*** *** +*** Wall time : 0h 0m 15s *** +*** *** +*************************************************************************** + + diff --git a/tests/li_scf_svwn5_libxc/test b/tests/li_scf_svwn5_libxc/test new file mode 100755 index 000000000..5aff74384 --- /dev/null +++ b/tests/li_scf_svwn5_libxc/test @@ -0,0 +1,26 @@ +#!/usr/bin/env python3 + +import sys +from pathlib import Path + +sys.path.append(str(Path(__file__).resolve().parents[1])) + +from tester import * # isort:skip + +options = script_cli() + +filters = { + SUM_OCCUPIED: rel_tolerance(1.0e-7), + E_KIN: rel_tolerance(1.0e-7), + E_EN: rel_tolerance(1.0e-7), + E_EE: rel_tolerance(1.0e-7), + E_X: rel_tolerance(1.0e-7), + E_XC: rel_tolerance(1.0e-7), + E_EEXT: rel_tolerance(1.0e-7), + E_NEXT: rel_tolerance(1.0e-7), + E_EL: rel_tolerance(1.0e-7), +} + +ierr = run(options, input_file="li", filters=filters) + +sys.exit(ierr)