1 ////////////////////////////////////////////////////////////////////////////////////// 2 // This file is distributed under the University of Illinois/NCSA Open Source License. 3 // See LICENSE file in top directory for details. 4 // 5 // Copyright (c) 2016 Jeongnim Kim and QMCPACK developers. 6 // 7 // File developed by: Mark Dewing, markdewing@gmail.com, University of Illinois at Urbana-Champaign 8 // 9 // File created by: Mark Dewing, markdewing@gmail.com, University of Illinois at Urbana-Champaign 10 ////////////////////////////////////////////////////////////////////////////////////// 11 12 13 #include "catch.hpp" 14 #include "Numerics/GaussianBasisSet.h" 15 16 #include <stdio.h> 17 #include <string> 18 19 using std::string; 20 21 typedef OHMMS_PRECISION real_type; 22 23 namespace qmcplusplus 24 { 25 TEST_CASE("Basic Gaussian", "[numerics]") 26 { 27 GaussianCombo<real_type>::BasicGaussian g1(3.0, 1.0); 28 29 real_type alpha = 3.0; 30 g1.reset(alpha, 1.0); 31 32 real_type r = 1.2; 33 34 real_type f = g1.f(r * r); 35 //std::cout << "f = " << f << std::endl << std::endl; 36 REQUIRE(f == Approx(0.0132998835424438)); 37 38 real_type df = g1.df(r, r * r); 39 //std::cout << "df = " << df << std::endl << std::endl; 40 REQUIRE(df == Approx(-0.0957591615055951)); 41 42 df = 0.0; 43 real_type ddf = 0.0; 44 f = g1.evaluate(r, r * r, df, ddf); 45 //std::cout << "f = " << f << " " << df << " " << ddf << std::endl; 46 REQUIRE(f == Approx(0.0132998835424438)); 47 REQUIRE(df == Approx(-0.0957591615055951)); 48 REQUIRE(ddf == Approx(0.609666661585622)); 49 50 df = 0.0; 51 ddf = 0.0; 52 real_type d3f = 0.0; 53 f = g1.evaluate(r, r * r, df, ddf, d3f); 54 //std::cout << "f = " << f << " " << df << " " << ddf << " " << d3f << std::endl; 55 REQUIRE(f == Approx(0.0132998835424438)); 56 REQUIRE(df == Approx(-0.0957591615055951)); 57 REQUIRE(ddf == Approx(0.609666661585622)); 58 REQUIRE(d3f == Approx(-3.24049002534934)); 59 } 60 61 TEST_CASE("Gaussian Combo", "[numerics]") 62 { 63 GaussianCombo<real_type> gc(0, false); 64 65 // STO-3G for H 66 gc.addGaussian(0.15432897, 3.42525091); 67 gc.addGaussian(0.53532814, 0.62391373); 68 gc.addGaussian(0.44463454, 0.16885540); 69 70 REQUIRE(gc.size() == 3); 71 72 real_type r = 1.3; 73 real_type f = gc.f(r); 74 //std::cout << "f = " << f << std::endl << std::endl; 75 REQUIRE(f == Approx(0.556240444149480)); 76 77 f = gc.evaluate(r, 1.0 / r); 78 REQUIRE(f == Approx(0.556240444149480)); 79 80 real_type df = gc.df(r); 81 REQUIRE(df == Approx(-0.661028435778766)); 82 83 gc.evaluateAll(r, 1.0 / r); 84 REQUIRE(gc.Y == Approx(0.556240444149480)); 85 REQUIRE(gc.dY == Approx(-0.661028435778766)); 86 REQUIRE(gc.d2Y == Approx(0.643259180749128)); 87 88 gc.evaluateWithThirdDeriv(r, 1.0 / r); 89 REQUIRE(gc.Y == Approx(0.556240444149480)); 90 REQUIRE(gc.dY == Approx(-0.661028435778766)); 91 REQUIRE(gc.d2Y == Approx(0.643259180749128)); 92 REQUIRE(gc.d3Y == Approx(-0.896186412781167)); 93 } 94 95 TEST_CASE("Gaussian Combo P", "[numerics]") 96 { 97 GaussianCombo<real_type> gc(1, false); 98 99 // cc-pVDZ for C, the 2P basis 100 gc.addGaussian(0.0381090, 9.4390000); 101 gc.addGaussian(0.2094800, 2.0020000); 102 gc.addGaussian(0.5085570, 0.5456000); 103 104 REQUIRE(gc.size() == 3); 105 106 real_type r = 1.3; 107 real_type f = gc.f(r); 108 109 REQUIRE(f == Approx(0.326057642350121)); 110 111 f = gc.evaluate(r, 1.0 / r); 112 REQUIRE(f == Approx(0.326057642350121)); 113 114 real_type df = gc.df(r); 115 REQUIRE(df == Approx(-0.649531407846947)); 116 117 gc.evaluateAll(r, 1.0 / r); 118 REQUIRE(gc.Y == Approx(0.326057642350121)); 119 REQUIRE(gc.dY == Approx(-0.649531407846947)); 120 REQUIRE(gc.d2Y == Approx(1.39522444199589)); 121 122 gc.evaluateWithThirdDeriv(r, 1.0 / r); 123 REQUIRE(gc.Y == Approx(0.326057642350121)); 124 REQUIRE(gc.dY == Approx(-0.649531407846947)); 125 REQUIRE(gc.d2Y == Approx(1.39522444199589)); 126 REQUIRE(gc.d3Y == Approx(-3.38467690038774)); 127 } 128 129 TEST_CASE("Gaussian Combo D", "[numerics]") 130 { 131 GaussianCombo<real_type> gc(2, false); 132 133 // cc-pVDZ for C, the 3D basis 134 gc.addGaussian(1.0, 0.5500000); 135 136 REQUIRE(gc.size() == 1); 137 138 real_type r = 1.3; 139 real_type f = gc.f(r); 140 141 REQUIRE(f == Approx(0.361815669819519)); 142 143 f = gc.evaluate(r, 1.0 / r); 144 REQUIRE(f == Approx(0.361815669819519)); 145 146 real_type df = gc.df(r); 147 REQUIRE(df == Approx(-0.517396407841913)); 148 149 gc.evaluateAll(r, 1.0 / r); 150 REQUIRE(gc.Y == Approx(0.361815669819519)); 151 REQUIRE(gc.dY == Approx(-0.517396407841913)); 152 REQUIRE(gc.d2Y == Approx(0.341879626412464)); 153 154 gc.evaluateWithThirdDeriv(r, 1.0 / r); 155 REQUIRE(gc.Y == Approx(0.361815669819519)); 156 REQUIRE(gc.dY == Approx(-0.517396407841913)); 157 REQUIRE(gc.d2Y == Approx(0.341879626412464)); 158 REQUIRE(gc.d3Y == Approx(0.649384231482385)); 159 } 160 } // namespace qmcplusplus 161