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: Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign
8 //                    Mark A. Berrill, berrillma@ornl.gov, Oak Ridge National Laboratory
9 //
10 // File created by: Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign
11 //////////////////////////////////////////////////////////////////////////////////////
12 
13 
14 #include "TwoDEwaldHandler.h"
15 
16 namespace qmcplusplus
17 {
initBreakup(ParticleSet & ref)18 void TwoDEwaldHandler::initBreakup(ParticleSet& ref)
19 {
20   LR_rc = ref.Lattice.LR_rc;
21   LR_kc = ref.Lattice.LR_kc;
22   Sigma = LR_kc;
23   //determine the sigma
24   while (erfc(Sigma) / LR_rc > 1e-16)
25   {
26     Sigma += 0.1;
27   }
28   app_log() << "   TwoDEwaldHandler Sigma/LR_rc = " << Sigma;
29   Sigma /= LR_rc;
30   app_log() << "  Sigma=" << Sigma << std::endl;
31   Volume = ref.Lattice.Volume;
32   fillFk(ref.SK->KLists);
33 }
34 
TwoDEwaldHandler(const TwoDEwaldHandler & aLR,ParticleSet & ref)35 TwoDEwaldHandler::TwoDEwaldHandler(const TwoDEwaldHandler& aLR, ParticleSet& ref)
36     : LRHandlerBase(aLR), Sigma(aLR.Sigma), Volume(aLR.Volume)
37 {}
38 
fillFk(KContainer & KList)39 void TwoDEwaldHandler::fillFk(KContainer& KList)
40 {
41   Fk.resize(KList.kpts_cart.size());
42   const std::vector<int>& kshell(KList.kshell);
43   MaxKshell = kshell.size() - 1;
44   Fk_symm.resize(MaxKshell);
45   kMag.resize(MaxKshell);
46   RealType knorm           = 2.0 * M_PI / Volume;
47   RealType oneovertwosigma = 1.0 / (2.0 * Sigma);
48   for (int ks = 0, ki = 0; ks < Fk_symm.size(); ks++)
49   {
50     kMag[ks]    = std::sqrt(KList.ksq[ki]);
51     RealType uk = knorm * erfc(kMag[ks] * oneovertwosigma) / kMag[ks];
52     Fk_symm[ks] = uk;
53     while (ki < KList.kshell[ks + 1] && ki < Fk.size())
54       Fk[ki++] = uk;
55     //       app_log()<<kMag[ks]<<" "<<uk<< std::endl;
56   }
57   app_log().flush();
58 }
59 } // namespace qmcplusplus
60