QMCPACK
test_LCAO_diamondC_2x1x1.cpp
Go to the documentation of this file.
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) 2023 QMCPACK developers.
6 //
7 // File developed by: Kevin Gasperich, kgasperich@anl.gov, Argonne National Laboratory
8 //
9 // File created by: Kevin Gasperich, kgasperich@anl.gov, Argonne National Laboratory
10 //////////////////////////////////////////////////////////////////////////////////////
11 
12 
13 #include "catch.hpp"
14 
15 #include "ParticleIO/LatticeIO.h"
16 #include "OhmmsData/Libxml2Doc.h"
17 #include "OhmmsPETE/OhmmsMatrix.h"
18 #include "Particle/ParticleSet.h"
19 #include "QMCWaveFunctions/WaveFunctionComponent.h"
20 #include "LCAO/LCAOrbitalBuilder.h"
21 #include <ResourceCollection.h>
22 #include "QMCHamiltonians/NLPPJob.h"
23 #include "DistanceTable.h"
24 
25 #include <stdio.h>
26 #include <string>
27 #include <limits>
28 
29 using std::string;
30 
31 namespace qmcplusplus
32 {
33 
34 template<typename DT>
35 using OffloadVector = Vector<DT, OffloadPinnedAllocator<DT>>;
36 
37 void test_LCAO_DiamondC_2x1x1_real(const bool useOffload)
38 {
39  using VT = SPOSet::ValueType;
40  Communicate* c = OHMMS::Controller;
41 
42  const char* particles = R"(<simulationcell>
43  <parameter name="lattice" units="bohr">
44  6.7463223 6.7463223 0.0
45  0.0 3.37316115 3.37316115
46  3.37316115 0.0 3.37316115
47  </parameter>
48  <parameter name="bconds">
49  p p p
50  </parameter>
51  </simulationcell>
52 )";
53 
54  Libxml2Document doc_lattice;
55  REQUIRE(doc_lattice.parseFromString(particles));
56 
57  // read lattice
58  ParticleSet::ParticleLayout lattice;
59  LatticeParser lp(lattice);
60  lp.put(doc_lattice.getRoot());
61  lattice.print(app_log(), 0);
62 
63  SimulationCell simcell(lattice);
64  ParticleSet ions_(simcell);
65 
66  ions_.setName("ion0");
67  ions_.create({4});
68  ions_.R[0] = {0.0, 0.0, 0.0};
69  ions_.R[1] = {1.686580575, 1.686580575, 1.686580575};
70  ions_.R[2] = {3.37316115, 3.37316115, 0.0};
71  ions_.R[3] = {5.059741726, 5.059741726, 1.686580575};
72  SpeciesSet& ispecies = ions_.getSpeciesSet();
73  const int Cidx = ispecies.addSpecies("C");
74 
75  ions_.print(app_log());
76  ions_.update(); // propagate SoA.
77 
78  ParticleSet elec_(simcell);
79  elec_.setName("elec");
80  elec_.create({8, 8});
81  elec_.R[0] = {0.0, 1.0, 0.0};
82  elec_.R[1] = {0.0, 1.1, 0.0};
83  elec_.R[2] = {0.0, 1.2, 0.0};
84  elec_.R[3] = {0.0, 1.3, 0.0};
85 
86  SpeciesSet& tspecies = elec_.getSpeciesSet();
87  const int upIdx = tspecies.addSpecies("u");
88  const int downIdx = tspecies.addSpecies("d");
89  const int chargeIdx = tspecies.addAttribute("charge");
90  tspecies(chargeIdx, upIdx) = -1;
91  tspecies(chargeIdx, downIdx) = -1;
92  const int ei_table_index = elec_.addTable(ions_);
93 
94  // diamondC_2x1x1
95  // from tests/solids/diamondC_2x1x1-Gaussian_pp/C_Diamond-Twist0.wfj.xml
96  const std::string wf_xml_str = R"(
97  <sposet_collection type="molecularorbital" name="LCAOBSet" source="ion0" transform="yes" twist="0 0 0" href="C_Diamond_2x1x1-Gaussian.h5" PBCimages="5 5 5" gpu="no">
98  <basisset name="LCAOBSet" key="GTO" transform="yes">
99  <grid type="log" ri="1.e-6" rf="1.e2" npts="1001"/>
100  </basisset>
101  <sposet name="spoud" size="8">
102  <occupation mode="ground"/>
103  <coefficient size="116" spindataset="0"/>
104  </sposet>
105  </sposet_collection>
106  )";
107  const std::string wf_omp_xml_str = R"(
108  <sposet_collection type="molecularorbital" name="LCAOBSet" source="ion0" transform="yes" twist="0 0 0" href="C_Diamond_2x1x1-Gaussian.h5" PBCimages="5 5 5" gpu="omptarget">
109  <basisset name="LCAOBSet" key="GTO" transform="yes">
110  <grid type="log" ri="1.e-6" rf="1.e2" npts="1001"/>
111  </basisset>
112  <sposet name="spoud" size="8">
113  <occupation mode="ground"/>
114  <coefficient size="116" spindataset="0"/>
115  </sposet>
116  </sposet_collection>
117  )";
118  Libxml2Document doc;
119  bool okay = doc.parseFromString(useOffload ? wf_omp_xml_str : wf_xml_str);
120  REQUIRE(okay);
121 
122  xmlNodePtr root = doc.getRoot();
123  xmlNodePtr bset_xml = xmlFirstElementChild(root);
124  xmlNodePtr sposet_xml = xmlNextElementSibling(bset_xml);
125 
126  LCAOrbitalBuilder lcaoSet(elec_, ions_, c, root);
127  auto spo = lcaoSet.createSPOSetFromXML(sposet_xml);
128  REQUIRE(spo);
129  auto& lcao_spos = dynamic_cast<const LCAOrbitalSet&>(*spo);
130  CHECK(!lcao_spos.isIdentity());
131 
132  const int norb = spo->getOrbitalSetSize();
133  app_log() << "norb: " << norb << std::endl;
134 
135  // test batched interfaces with 2 walkers
136  const size_t nw = 2;
137 
138  ParticleSet elec_2(elec_);
139  // interchange positions
140  elec_2.R[0] = elec_.R[1];
141  elec_2.R[1] = elec_.R[0];
142 
143  RefVectorWithLeader<ParticleSet> p_list(elec_, {elec_, elec_2});
144 
145  std::unique_ptr<SPOSet> spo_2(spo->makeClone());
146  RefVectorWithLeader<SPOSet> spo_list(*spo, {*spo, *spo_2});
147 
148  ResourceCollection pset_res("test_pset_res");
149  ResourceCollection spo_res("test_spo_res");
150 
151  elec_.createResource(pset_res);
152  spo->createResource(spo_res);
153 
154  ResourceCollectionTeamLock<ParticleSet> mw_pset_lock(pset_res, p_list);
155  ResourceCollectionTeamLock<SPOSet> mw_sposet_lock(spo_res, spo_list);
156  ParticleSet::mw_update(p_list);
157 
158  SECTION("LCAOrbitalSet::mw_evaluateVGL")
159  {
160  SPOSet::ValueVector psiref_0(norb);
161  SPOSet::GradVector dpsiref_0(norb);
162  SPOSet::ValueVector d2psiref_0(norb);
163  SPOSet::ValueVector psiref_1(norb);
164  SPOSet::GradVector dpsiref_1(norb);
165  SPOSet::ValueVector d2psiref_1(norb);
166 
167  spo_2->evaluateVGL(elec_2, 0, psiref_1, dpsiref_1, d2psiref_1);
168  spo->evaluateVGL(elec_, 0, psiref_0, dpsiref_0, d2psiref_0);
169  SECTION("single-walker VGL")
170  {
171  CHECK(Approx(std::real(psiref_0[0])) == 0.10203360788082);
172  CHECK(Approx(std::real(d2psiref_0[0])) == -0.14269632514649);
173  CHECK(Approx(std::real(dpsiref_0[0][0])) == 0.00031796000022);
174  CHECK(Approx(std::real(dpsiref_0[0][1])) == -0.01951178212495);
175  CHECK(Approx(std::real(dpsiref_0[0][2])) == -0.00031738324507);
176  CHECK(Approx(std::real(psiref_0[1])) == 0.14111282090404);
177  CHECK(Approx(std::real(d2psiref_0[1])) == -0.28230884342631);
178  CHECK(Approx(std::real(dpsiref_0[1][0])) == 0.01205368790709);
179  CHECK(Approx(std::real(dpsiref_0[1][1])) == -0.04876640907938);
180  CHECK(Approx(std::real(dpsiref_0[1][2])) == -0.01205402562448);
181 
182  CHECK(Approx(std::real(psiref_1[0])) == 0.09954792826469);
183  CHECK(Approx(std::real(d2psiref_1[0])) == -0.10799468581785);
184  CHECK(Approx(std::real(dpsiref_1[0][0])) == 0.00024577684629);
185  CHECK(Approx(std::real(dpsiref_1[0][1])) == -0.02943596305516);
186  CHECK(Approx(std::real(dpsiref_1[0][2])) == -0.00024512723822);
187  CHECK(Approx(std::real(psiref_1[1])) == 0.13558495733438);
188  CHECK(Approx(std::real(d2psiref_1[1])) == -0.21885179829267);
189  CHECK(Approx(std::real(dpsiref_1[1][0])) == 0.00738771300147);
190  CHECK(Approx(std::real(dpsiref_1[1][1])) == -0.06080141726019);
191  CHECK(Approx(std::real(dpsiref_1[1][2])) == -0.00738811343748);
192  }
193 
194  SPOSet::ValueVector psi_v_1(norb);
195  SPOSet::ValueVector psi_v_2(norb);
196  RefVector<SPOSet::ValueVector> psi_v_list{psi_v_1, psi_v_2};
197  spo->mw_evaluateValue(spo_list, p_list, 0, psi_v_list);
198  SECTION("multi-walker V")
199  {
200  CHECK(Approx(std::real(psi_v_list[0].get()[0])) == 0.10203360788082);
201  CHECK(Approx(std::real(psi_v_list[0].get()[1])) == 0.14111282090404);
202  CHECK(Approx(std::real(psi_v_list[1].get()[0])) == 0.09954792826469);
203  CHECK(Approx(std::real(psi_v_list[1].get()[1])) == 0.13558495733438);
204  }
205 
206 
207  SPOSet::ValueVector psi_1(norb);
208  SPOSet::GradVector dpsi_1(norb);
209  SPOSet::ValueVector d2psi_1(norb);
210  SPOSet::ValueVector psi_2(norb);
211  SPOSet::GradVector dpsi_2(norb);
212  SPOSet::ValueVector d2psi_2(norb);
213  RefVector<SPOSet::ValueVector> psi_list = {psi_1, psi_2};
214  RefVector<SPOSet::GradVector> dpsi_list = {dpsi_1, dpsi_2};
215  RefVector<SPOSet::ValueVector> d2psi_list = {d2psi_1, d2psi_2};
216  spo->mw_evaluateVGL(spo_list, p_list, 0, psi_list, dpsi_list, d2psi_list);
217  SECTION("multi-walker VGL")
218  {
219  CHECK(Approx(std::real(psi_list[0].get()[0])) == 0.10203360788082);
220  CHECK(Approx(std::real(d2psi_list[0].get()[0])) == -0.14269632514649);
221  CHECK(Approx(std::real(dpsi_list[0].get()[0][0])) == 0.00031796000022);
222  CHECK(Approx(std::real(dpsi_list[0].get()[0][1])) == -0.01951178212495);
223  CHECK(Approx(std::real(dpsi_list[0].get()[0][2])) == -0.00031738324507);
224  CHECK(Approx(std::real(psi_list[0].get()[1])) == 0.14111282090404);
225  CHECK(Approx(std::real(d2psi_list[0].get()[1])) == -0.28230884342631);
226  CHECK(Approx(std::real(dpsi_list[0].get()[1][0])) == 0.01205368790709);
227  CHECK(Approx(std::real(dpsi_list[0].get()[1][1])) == -0.04876640907938);
228  CHECK(Approx(std::real(dpsi_list[0].get()[1][2])) == -0.01205402562448);
229 
230  CHECK(Approx(std::real(psi_list[1].get()[0])) == 0.09954792826469);
231  CHECK(Approx(std::real(d2psi_list[1].get()[0])) == -0.10799468581785);
232  CHECK(Approx(std::real(dpsi_list[1].get()[0][0])) == 0.00024577684629);
233  CHECK(Approx(std::real(dpsi_list[1].get()[0][1])) == -0.02943596305516);
234  CHECK(Approx(std::real(dpsi_list[1].get()[0][2])) == -0.00024512723822);
235  CHECK(Approx(std::real(psi_list[1].get()[1])) == 0.13558495733438);
236  CHECK(Approx(std::real(d2psi_list[1].get()[1])) == -0.21885179829267);
237  CHECK(Approx(std::real(dpsi_list[1].get()[1][0])) == 0.00738771300147);
238  CHECK(Approx(std::real(dpsi_list[1].get()[1][1])) == -0.06080141726019);
239  CHECK(Approx(std::real(dpsi_list[1].get()[1][2])) == -0.00738811343748);
240  }
241 
242 
243  SECTION("compare MW/SW V/VGL")
244  {
245  for (size_t iorb = 0; iorb < norb; iorb++)
246  {
247  CHECK(Approx(std::real(psi_list[0].get()[iorb])) == std::real(psiref_0[iorb]));
248  CHECK(Approx(std::real(psi_list[1].get()[iorb])) == std::real(psiref_1[iorb]));
249  CHECK(Approx(std::real(d2psi_list[0].get()[iorb])) == std::real(d2psiref_0[iorb]));
250  CHECK(Approx(std::real(d2psi_list[1].get()[iorb])) == std::real(d2psiref_1[iorb]));
251  for (size_t idim = 0; idim < SPOSet::DIM; idim++)
252  {
253  CHECK(Approx(std::real(dpsi_list[0].get()[iorb][idim])) == std::real(dpsiref_0[iorb][idim]));
254  CHECK(Approx(std::real(dpsi_list[1].get()[iorb][idim])) == std::real(dpsiref_1[iorb][idim]));
255  }
256  }
257  for (size_t iorb = 0; iorb < norb; iorb++)
258  for (size_t iw = 0; iw < nw; iw++)
259  CHECK(Approx(std::real(psi_v_list[iw].get()[iorb])) == std::real(psi_list[iw].get()[iorb]));
260 #ifdef QMC_COMPLEX
261  for (size_t iorb = 0; iorb < norb; iorb++)
262  {
263  CHECK(Approx(std::imag(psi_list[0].get()[iorb])) == std::imag(psiref_0[iorb]));
264  CHECK(Approx(std::imag(psi_list[1].get()[iorb])) == std::imag(psiref_1[iorb]));
265  CHECK(Approx(std::imag(d2psi_list[0].get()[iorb])) == std::imag(d2psiref_0[iorb]));
266  CHECK(Approx(std::imag(d2psi_list[1].get()[iorb])) == std::imag(d2psiref_1[iorb]));
267  for (size_t idim = 0; idim < SPOSet::DIM; idim++)
268  {
269  CHECK(Approx(std::imag(dpsi_list[0].get()[iorb][idim])) == std::imag(dpsiref_0[iorb][idim]));
270  CHECK(Approx(std::imag(dpsi_list[1].get()[iorb][idim])) == std::imag(dpsiref_1[iorb][idim]));
271  }
272  }
273  for (size_t iorb = 0; iorb < norb; iorb++)
274  for (size_t iw = 0; iw < nw; iw++)
275  CHECK(Approx(std::imag(psi_v_list[iw].get()[iorb])) == std::imag(psi_list[iw].get()[iorb]));
276 #endif
277  }
278  // app_log() << "vgl_refvalues: \n" << std::setprecision(14);
279  // for (int iorb = 0; iorb < 2; iorb++)
280  // {
281  // app_log() << "CHECK(Approx(std::real(psiref_0[" << iorb << "])) == " << psiref_0[iorb] << ");\n";
282  // app_log() << "CHECK(Approx(std::real(d2psiref_0[" << iorb << "])) == " << d2psiref_0[iorb] << ");\n";
283  // for (int idim = 0; idim < 3; idim++)
284  // app_log() << "CHECK(Approx(std::real(dpsiref_0[" << iorb << "][" << idim << "])) == " << dpsiref_0[iorb][idim] << ");\n";
285  // }
286  // for (int iorb = 0; iorb < 2; iorb++)
287  // {
288  // app_log() << "CHECK(Approx(std::real(psiref_1[" << iorb << "])) == " << psiref_1[iorb] << ");\n";
289  // app_log() << "CHECK(Approx(std::real(d2psiref_1[" << iorb << "])) == " << d2psiref_1[iorb] << "));\n";
290  // for (int idim = 0; idim < 3; idim++)
291  // app_log() << "CHECK(Approx(std::real(dpsiref_1[" << iorb << "][" << idim << "])) == " << dpsiref_1[iorb][idim] << ");\n";
292  // }
293 
294  // fill invrow with dummy data for each walker
295  OffloadVector<SPOSet::ValueType> psiMinv_data_0(norb), psiMinv_data_1(norb);
296  LCAOrbitalSet::ValueVector psiMinv_ref_0(psiMinv_data_0.data(), norb), psiMinv_ref_1(psiMinv_data_1.data(), norb);
297  for (int i = 0; i < norb; i++)
298  {
299  psiMinv_data_0[i] = 0.1 * i;
300  psiMinv_data_1[i] = 0.2 * i;
301  }
302  psiMinv_data_0.updateTo();
303  psiMinv_data_1.updateTo();
304 
305  std::vector<const SPOSet::ValueType*> invRow_ptr_list;
306  if (spo->isOMPoffload())
307  invRow_ptr_list = {psiMinv_data_0.device_data(), psiMinv_data_1.device_data()};
308  else
309  invRow_ptr_list = {psiMinv_data_0.data(), psiMinv_data_1.data()};
310 
311  SPOSet::OffloadMWVGLArray phi_vgl_v;
312  std::vector<SPOSet::ValueType> ratios(nw);
313  std::vector<SPOSet::GradType> grads(nw);
314  phi_vgl_v.resize(QMCTraits::DIM_VGL, nw, norb);
315 
316  spo->mw_evaluateVGLandDetRatioGrads(spo_list, p_list, 0, invRow_ptr_list, phi_vgl_v, ratios, grads);
317 
318  SECTION("multi-walker VGLandDetRatioGrads")
319  {
320  CHECK(Approx(std::real(ratios[0])) == 0.193096895);
321  CHECK(Approx(std::real(grads[0][0])) == 0.5969699486);
322  CHECK(Approx(std::real(grads[0][1])) == -0.776008772);
323  CHECK(Approx(std::real(grads[0][2])) == 0.7457653703);
324  CHECK(Approx(std::real(ratios[1])) == 0.3558058932);
325  CHECK(Approx(std::real(grads[1][0])) == 0.5857037763);
326  CHECK(Approx(std::real(grads[1][1])) == -0.8617132632);
327  CHECK(Approx(std::real(grads[1][2])) == 0.7299292825);
328  }
329  }
330 
331  SECTION("LCAOrbitalSet::mw_evaluateDetRatios")
332  {
333  // make VPs
334  const size_t nvp_ = 4;
335  const size_t nvp_2 = 3;
336  const std::vector<size_t> nvp_list = {nvp_, nvp_2};
337  VirtualParticleSet VP_(elec_, nvp_);
338  VirtualParticleSet VP_2(elec_2, nvp_2);
339 
340  // move VPs
341  std::vector<ParticleSet::SingleParticlePos> newpos_vp_(nvp_);
342  std::vector<ParticleSet::SingleParticlePos> newpos_vp_2(nvp_2);
343  for (int i = 0; i < nvp_; i++)
344  {
345  newpos_vp_[i][0] = 0.1 * i;
346  newpos_vp_[i][1] = 0.2 * i;
347  newpos_vp_[i][2] = 0.3 * i;
348  }
349  for (int i = 0; i < nvp_2; i++)
350  {
351  newpos_vp_2[i][0] = 0.2 * i;
352  newpos_vp_2[i][1] = 0.3 * i;
353  newpos_vp_2[i][2] = 0.4 * i;
354  }
355 
356  const auto& ei_table_ = elec_.getDistTableAB(ei_table_index);
357  const auto& ei_table_2 = elec_2.getDistTableAB(ei_table_index);
358  // make virtual move of elec 0, reference ion 1
359  NLPPJob<SPOSet::RealType> job_(1, 0, ei_table_.getDistances()[0][1], -ei_table_.getDisplacements()[0][1]);
360  // make virtual move of elec 1, reference ion 3
361  NLPPJob<SPOSet::RealType> job_2(3, 1, ei_table_2.getDistances()[1][3], -ei_table_2.getDisplacements()[1][3]);
362 
363  // make VP refvec
364  RefVectorWithLeader<VirtualParticleSet> vp_list(VP_, {VP_, VP_2});
365 
366  ResourceCollection vp_res("test_vp_res");
367  VP_.createResource(vp_res);
368  ResourceCollectionTeamLock<VirtualParticleSet> mw_vpset_lock(vp_res, vp_list);
369  VirtualParticleSet::mw_makeMoves(vp_list, p_list, {newpos_vp_, newpos_vp_2}, {job_, job_2}, false);
370 
371  // fill invrow with dummy data for each walker
372  OffloadVector<SPOSet::ValueType> psiMinv_data_0(norb), psiMinv_data_1(norb);
373  LCAOrbitalSet::ValueVector psiMinv_ref_0(psiMinv_data_0.data(), norb), psiMinv_ref_1(psiMinv_data_1.data(), norb);
374  for (int i = 0; i < norb; i++)
375  {
376  psiMinv_data_0[i] = 0.1 * i;
377  psiMinv_data_1[i] = 0.2 * i;
378  }
379  psiMinv_data_0.updateTo();
380  psiMinv_data_1.updateTo();
381 
382  std::vector<const SPOSet::ValueType*> invRow_ptr_list;
383  if (spo->isOMPoffload())
384  invRow_ptr_list = {psiMinv_data_0.device_data(), psiMinv_data_1.device_data()};
385  else
386  invRow_ptr_list = {psiMinv_data_0.data(), psiMinv_data_1.data()};
387 
388  // ratios_list
389  std::vector<std::vector<SPOSet::ValueType>> ratios_list(nw);
390  for (size_t iw = 0; iw < nw; iw++)
391  ratios_list[iw].resize(nvp_list[iw]);
392 
393  // just need dummy refvec with correct size
394  SPOSet::ValueVector tmp_psi0(norb), tmp_psi1(norb);
395  RefVector<SPOSet::ValueVector> tmp_psi_list{tmp_psi0, tmp_psi1};
396  spo->mw_evaluateDetRatios(spo_list, RefVectorWithLeader<const VirtualParticleSet>(VP_, {VP_, VP_2}), tmp_psi_list,
397  invRow_ptr_list, ratios_list);
398 
399  std::vector<SPOSet::ValueType> ratios_ref_0(nvp_);
400  std::vector<SPOSet::ValueType> ratios_ref_1(nvp_2);
401  // single-walker functions for reference
402  spo->evaluateDetRatios(VP_, tmp_psi0, psiMinv_ref_0, ratios_ref_0);
403  spo_2->evaluateDetRatios(VP_2, tmp_psi1, psiMinv_ref_1, ratios_ref_1);
404 
405  for (int ivp = 0; ivp < nvp_; ivp++)
406  CHECK(Approx(std::real(ratios_list[0][ivp])) == std::real(ratios_ref_0[ivp]));
407  for (int ivp = 0; ivp < nvp_2; ivp++)
408  CHECK(Approx(std::real(ratios_list[1][ivp])) == std::real(ratios_ref_1[ivp]));
409 #ifdef QMC_COMPLEX
410  for (int ivp = 0; ivp < nvp_; ivp++)
411  {
412  CHECK(Approx(std::imag(ratios_list[0][ivp])) == std::imag(ratios_ref_0[ivp]));
413  CHECK(Approx(std::imag(ratios_list[0][ivp])) == 0.0);
414  }
415  for (int ivp = 0; ivp < nvp_2; ivp++)
416  {
417  CHECK(Approx(std::imag(ratios_list[1][ivp])) == std::imag(ratios_ref_1[ivp]));
418  CHECK(Approx(std::imag(ratios_list[1][ivp])) == 0.0);
419  }
420 #endif
421 
422  CHECK(Approx(std::real(ratios_list[0][0])) == 0.19309684969511);
423  CHECK(Approx(std::real(ratios_list[0][1])) == 0.19743141486366);
424  CHECK(Approx(std::real(ratios_list[0][2])) == 0.17884881050205);
425  CHECK(Approx(std::real(ratios_list[0][3])) == 0.15105783567230);
426  CHECK(Approx(std::real(ratios_list[1][0])) == 0.38619369939021);
427  CHECK(Approx(std::real(ratios_list[1][1])) == 0.38429955941922);
428  CHECK(Approx(std::real(ratios_list[1][2])) == 0.32071997896196);
429 
430  CHECK(Approx(std::real(ratios_ref_0[0])) == 0.1930968497);
431  CHECK(Approx(std::real(ratios_ref_0[1])) == 0.1974314149);
432  CHECK(Approx(std::real(ratios_ref_0[2])) == 0.1788488105);
433  CHECK(Approx(std::real(ratios_ref_0[3])) == 0.1510578357);
434  CHECK(Approx(std::real(ratios_ref_1[0])) == 0.3861936994);
435  CHECK(Approx(std::real(ratios_ref_1[1])) == 0.3842995594);
436  CHECK(Approx(std::real(ratios_ref_1[2])) == 0.3207199790);
437 
438  //// print SW ref values
439  // for (int ivp = 0; ivp < nvp_; ivp++)
440  // app_log() << "CHECK(Approx(std::real(ratios_ref_0[" << ivp << "])) == " << std::real(ratios_ref_0[ivp]) << ");\n";
441  // for (int ivp = 0; ivp < nvp_2; ivp++)
442  // app_log() << "CHECK(Approx(std::real(ratios_ref_1[" << ivp << "])) == " << std::real(ratios_ref_1[ivp]) << ");\n";
443 
444  //// print MW ref values
445  // for (int iw = 0; iw < nw; iw++)
446  // for (int ivp = 0; ivp < nvp_list[iw]; ivp++)
447  // app_log() << "CHECK(Approx(std::real(ratios_list[" << iw << "][" << ivp << "])) == " << std::real(ratios_list[iw][ivp]) << ");\n";
448  // for (int iw = 0; iw < nw; iw++)
449  // for (int ivp = 0; ivp < nvp_list[iw]; ivp++)
450  // app_log() << "CHECK(Approx(std::imag(ratios_list[" << iw << "][" << ivp << "])) == " << std::imag(ratios_list[iw][ivp]) << ");\n";
451  }
452 }
453 
454 void test_LCAO_DiamondC_2x1x1_cplx(const bool useOffload)
455 {
456  using VT = SPOSet::ValueType;
457  Communicate* c = OHMMS::Controller;
458 
459  const char* particles = R"(<simulationcell>
460  <parameter name="lattice" units="bohr">
461  6.7463223 6.7463223 0.0
462  0.0 3.37316115 3.37316115
463  3.37316115 0.0 3.37316115
464  </parameter>
465  <parameter name="bconds">
466  p p p
467  </parameter>
468  </simulationcell>
469 )";
470 
471  Libxml2Document doc_lattice;
472  REQUIRE(doc_lattice.parseFromString(particles));
473 
474  // read lattice
475  ParticleSet::ParticleLayout lattice;
476  LatticeParser lp(lattice);
477  lp.put(doc_lattice.getRoot());
478  lattice.print(app_log(), 0);
479 
480  SimulationCell simcell(lattice);
481  ParticleSet ions_(simcell);
482 
483  ions_.setName("ion0");
484  ions_.create({4});
485  ions_.R[0] = {0.0, 0.0, 0.0};
486  ions_.R[1] = {1.686580575, 1.686580575, 1.686580575};
487  ions_.R[2] = {3.37316115, 3.37316115, 0.0};
488  ions_.R[3] = {5.059741726, 5.059741726, 1.686580575};
489  SpeciesSet& ispecies = ions_.getSpeciesSet();
490  const int Cidx = ispecies.addSpecies("C");
491 
492  ions_.print(app_log());
493  ions_.update(); // propagate SoA.
494 
495  ParticleSet elec_(simcell);
496  elec_.setName("elec");
497  elec_.create({8, 8});
498  elec_.R[0] = {0.0, 1.0, 0.0};
499  elec_.R[1] = {0.0, 1.1, 0.0};
500  elec_.R[2] = {0.0, 1.2, 0.0};
501  elec_.R[3] = {0.0, 1.3, 0.0};
502 
503  SpeciesSet& tspecies = elec_.getSpeciesSet();
504  const int upIdx = tspecies.addSpecies("u");
505  const int downIdx = tspecies.addSpecies("d");
506  const int chargeIdx = tspecies.addAttribute("charge");
507  tspecies(chargeIdx, upIdx) = -1;
508  tspecies(chargeIdx, downIdx) = -1;
509  const int ei_table_index = elec_.addTable(ions_);
510 
511  // diamondC_2x1x1
512  // from tests/solids/diamondC_2x1x1-Gaussian_pp_cplx/C_Diamond-tiled-cplx.wfj.xml
513  const std::string wf_xml_str = R"(
514  <sposet_collection type="molecularorbital" name="LCAOBSet" source="ion0" transform="yes" twist="0.07761248 0.07761248 -0.07761248" href="C_Diamond_2x1x1-Gaussian-tiled-cplx.h5" PBCimages="5 5 5" gpu="no">
515  <basisset name="LCAOBSet" key="GTO" transform="yes">
516  <grid type="log" ri="1.e-6" rf="1.e2" npts="1001"/>
517  </basisset>
518  <sposet name="spoud" size="8" >
519  <occupation mode="ground"/>
520  <coefficient size="52" spindataset="0"/>
521  </sposet>
522  </sposet_collection>
523  )";
524  const std::string wf_omp_xml_str = R"(
525  <sposet_collection type="molecularorbital" name="LCAOBSet" source="ion0" transform="yes" twist="0.07761248 0.07761248 -0.07761248" href="C_Diamond_2x1x1-Gaussian-tiled-cplx.h5" PBCimages="5 5 5" gpu="omptarget">
526  <basisset name="LCAOBSet" key="GTO" transform="yes">
527  <grid type="log" ri="1.e-6" rf="1.e2" npts="1001"/>
528  </basisset>
529  <sposet name="spoud" size="8" >
530  <occupation mode="ground"/>
531  <coefficient size="52" spindataset="0"/>
532  </sposet>
533  </sposet_collection>
534  )";
535  Libxml2Document doc;
536  bool okay = doc.parseFromString(useOffload ? wf_omp_xml_str : wf_xml_str);
537  REQUIRE(okay);
538 
539  xmlNodePtr root = doc.getRoot();
540  xmlNodePtr bset_xml = xmlFirstElementChild(root);
541  xmlNodePtr sposet_xml = xmlNextElementSibling(bset_xml);
542 
543  LCAOrbitalBuilder lcaoSet(elec_, ions_, c, root);
544  auto spo = lcaoSet.createSPOSetFromXML(sposet_xml);
545  REQUIRE(spo);
546  auto& lcao_spos = dynamic_cast<const LCAOrbitalSet&>(*spo);
547  CHECK(!lcao_spos.isIdentity());
548 
549  const int norb = spo->getOrbitalSetSize();
550  app_log() << "norb: " << norb << std::endl;
551 
552  // test batched interfaces with 2 walkers
553  const size_t nw = 2;
554 
555  ParticleSet elec_2(elec_);
556  // interchange positions
557  elec_2.R[0] = elec_.R[1];
558  elec_2.R[1] = elec_.R[0];
559 
560  RefVectorWithLeader<ParticleSet> p_list(elec_, {elec_, elec_2});
561 
562  std::unique_ptr<SPOSet> spo_2(spo->makeClone());
563  RefVectorWithLeader<SPOSet> spo_list(*spo, {*spo, *spo_2});
564 
565  ResourceCollection pset_res("test_pset_res");
566  ResourceCollection spo_res("test_spo_res");
567 
568  elec_.createResource(pset_res);
569  spo->createResource(spo_res);
570 
571  ResourceCollectionTeamLock<ParticleSet> mw_pset_lock(pset_res, p_list);
572  ResourceCollectionTeamLock<SPOSet> mw_sposet_lock(spo_res, spo_list);
573  ParticleSet::mw_update(p_list);
574 
575  SECTION("LCAOrbitalSet::mw_evaluateVGL")
576  {
577  SPOSet::ValueVector psiref_0(norb);
578  SPOSet::GradVector dpsiref_0(norb);
579  SPOSet::ValueVector d2psiref_0(norb);
580  SPOSet::ValueVector psiref_1(norb);
581  SPOSet::GradVector dpsiref_1(norb);
582  SPOSet::ValueVector d2psiref_1(norb);
583 
584  spo_2->evaluateVGL(elec_2, 0, psiref_1, dpsiref_1, d2psiref_1);
585  spo->evaluateVGL(elec_, 0, psiref_0, dpsiref_0, d2psiref_0);
586  SECTION("single-walker VGL")
587  {
588  CHECK(Approx(std::real(psiref_0[0])) == 0.10394953298732);
589  CHECK(Approx(std::real(d2psiref_0[0])) == -0.15355833615767);
590  CHECK(Approx(std::real(dpsiref_0[0][0])) == -0.00066360565262);
591  CHECK(Approx(std::real(dpsiref_0[0][1])) == -0.03040447635885);
592  CHECK(Approx(std::real(dpsiref_0[0][2])) == 0.00066400579544);
593  CHECK(Approx(std::real(psiref_0[1])) == 0.04007076321982);
594  CHECK(Approx(std::real(d2psiref_0[1])) == -0.13347762878346);
595  CHECK(Approx(std::real(dpsiref_0[1][0])) == 0.06427121424074);
596  CHECK(Approx(std::real(dpsiref_0[1][1])) == -0.03145499248870);
597  CHECK(Approx(std::real(dpsiref_0[1][2])) == 0.09703977123104);
598 
599  CHECK(Approx(std::real(psiref_1[0])) == 0.10041826765555);
600  CHECK(Approx(std::real(d2psiref_1[0])) == -0.11271657756448);
601  CHECK(Approx(std::real(dpsiref_1[0][0])) == -0.00056949645372);
602  CHECK(Approx(std::real(dpsiref_1[0][1])) == -0.03940982127946);
603  CHECK(Approx(std::real(dpsiref_1[0][2])) == 0.00056995118020);
604  CHECK(Approx(std::real(psiref_1[1])) == 0.03692120161253);
605  CHECK(Approx(std::real(d2psiref_1[1])) == -0.10225156633271);
606  CHECK(Approx(std::real(dpsiref_1[1][0])) == 0.05862729651642);
607  CHECK(Approx(std::real(dpsiref_1[1][1])) == -0.03141220770622);
608  CHECK(Approx(std::real(dpsiref_1[1][2])) == 0.08454924955444);
609 
610  CHECK(Approx(std::imag(psiref_0[0])) == 0.00920567821605);
611  CHECK(Approx(std::imag(d2psiref_0[0])) == -0.01897980702117);
612  CHECK(Approx(std::imag(dpsiref_0[0][0])) == 0.00599757336742);
613  CHECK(Approx(std::imag(dpsiref_0[0][1])) == 0.00009471083794);
614  CHECK(Approx(std::imag(dpsiref_0[0][2])) == -0.00599734141981);
615  CHECK(Approx(std::imag(psiref_0[1])) == -0.07297020296665);
616  CHECK(Approx(std::imag(d2psiref_0[1])) == 0.24972839540365);
617  CHECK(Approx(std::imag(dpsiref_0[1][0])) == -0.14785027769482);
618  CHECK(Approx(std::imag(dpsiref_0[1][1])) == 0.05546078377043);
619  CHECK(Approx(std::imag(dpsiref_0[1][2])) == -0.19276310593334);
620 
621  CHECK(Approx(std::imag(psiref_1[0])) == 0.00921878891696);
622  CHECK(Approx(std::imag(d2psiref_1[0])) == -0.01501853597019);
623  CHECK(Approx(std::imag(dpsiref_1[0][0])) == 0.00531194469490);
624  CHECK(Approx(std::imag(dpsiref_1[0][1])) == 0.00017539927704);
625  CHECK(Approx(std::imag(dpsiref_1[0][2])) == -0.00531168663225);
626  CHECK(Approx(std::imag(psiref_1[1])) == -0.06736883166719);
627  CHECK(Approx(std::imag(d2psiref_1[1])) == 0.19452200419216);
628  CHECK(Approx(std::imag(dpsiref_1[1][0])) == -0.13949226666533);
629  CHECK(Approx(std::imag(dpsiref_1[1][1])) == 0.05635953434634);
630  CHECK(Approx(std::imag(dpsiref_1[1][2])) == -0.17227552311686);
631  }
632 
633  SPOSet::ValueVector psi_v_1(norb);
634  SPOSet::ValueVector psi_v_2(norb);
635  RefVector<SPOSet::ValueVector> psi_v_list{psi_v_1, psi_v_2};
636  spo->mw_evaluateValue(spo_list, p_list, 0, psi_v_list);
637  SECTION("multi-walker V")
638  {
639  CHECK(Approx(std::real(psi_v_list[0].get()[0])) == 0.10394953298732);
640  CHECK(Approx(std::real(psi_v_list[0].get()[1])) == 0.04007076321982);
641  CHECK(Approx(std::real(psi_v_list[1].get()[0])) == 0.10041826765555);
642  CHECK(Approx(std::real(psi_v_list[1].get()[1])) == 0.03692120161253);
643 
644  CHECK(Approx(std::imag(psi_v_list[0].get()[0])) == 0.00920567821605);
645  CHECK(Approx(std::imag(psi_v_list[0].get()[1])) == -0.07297020296665);
646  CHECK(Approx(std::imag(psi_v_list[1].get()[0])) == 0.00921878891696);
647  CHECK(Approx(std::imag(psi_v_list[1].get()[1])) == -0.06736883166719);
648  }
649 
650 
651  SPOSet::ValueVector psi_1(norb);
652  SPOSet::GradVector dpsi_1(norb);
653  SPOSet::ValueVector d2psi_1(norb);
654  SPOSet::ValueVector psi_2(norb);
655  SPOSet::GradVector dpsi_2(norb);
656  SPOSet::ValueVector d2psi_2(norb);
657  RefVector<SPOSet::ValueVector> psi_list = {psi_1, psi_2};
658  RefVector<SPOSet::GradVector> dpsi_list = {dpsi_1, dpsi_2};
659  RefVector<SPOSet::ValueVector> d2psi_list = {d2psi_1, d2psi_2};
660  spo->mw_evaluateVGL(spo_list, p_list, 0, psi_list, dpsi_list, d2psi_list);
661  SECTION("multi-walker VGL")
662  {
663  CHECK(Approx(std::real(psi_list[0].get()[0])) == 0.10394953298732);
664  CHECK(Approx(std::real(d2psi_list[0].get()[0])) == -0.15355833615767);
665  CHECK(Approx(std::real(dpsi_list[0].get()[0][0])) == -0.00066360565262);
666  CHECK(Approx(std::real(dpsi_list[0].get()[0][1])) == -0.03040447635885);
667  CHECK(Approx(std::real(dpsi_list[0].get()[0][2])) == 0.00066400579544);
668  CHECK(Approx(std::real(psi_list[0].get()[1])) == 0.04007076321982);
669  CHECK(Approx(std::real(d2psi_list[0].get()[1])) == -0.13347762878346);
670  CHECK(Approx(std::real(dpsi_list[0].get()[1][0])) == 0.06427121424074);
671  CHECK(Approx(std::real(dpsi_list[0].get()[1][1])) == -0.03145499248870);
672  CHECK(Approx(std::real(dpsi_list[0].get()[1][2])) == 0.09703977123104);
673 
674  CHECK(Approx(std::real(psi_list[1].get()[0])) == 0.10041826765555);
675  CHECK(Approx(std::real(d2psi_list[1].get()[0])) == -0.11271657756448);
676  CHECK(Approx(std::real(dpsi_list[1].get()[0][0])) == -0.00056949645372);
677  CHECK(Approx(std::real(dpsi_list[1].get()[0][1])) == -0.03940982127946);
678  CHECK(Approx(std::real(dpsi_list[1].get()[0][2])) == 0.00056995118020);
679  CHECK(Approx(std::real(psi_list[1].get()[1])) == 0.03692120161253);
680  CHECK(Approx(std::real(d2psi_list[1].get()[1])) == -0.10225156633271);
681  CHECK(Approx(std::real(dpsi_list[1].get()[1][0])) == 0.05862729651642);
682  CHECK(Approx(std::real(dpsi_list[1].get()[1][1])) == -0.03141220770622);
683  CHECK(Approx(std::real(dpsi_list[1].get()[1][2])) == 0.08454924955444);
684 
685  CHECK(Approx(std::imag(psi_list[0].get()[0])) == 0.00920567821605);
686  CHECK(Approx(std::imag(d2psi_list[0].get()[0])) == -0.01897980702117);
687  CHECK(Approx(std::imag(dpsi_list[0].get()[0][0])) == 0.00599757336742);
688  CHECK(Approx(std::imag(dpsi_list[0].get()[0][1])) == 0.00009471083794);
689  CHECK(Approx(std::imag(dpsi_list[0].get()[0][2])) == -0.00599734141981);
690  CHECK(Approx(std::imag(psi_list[0].get()[1])) == -0.07297020296665);
691  CHECK(Approx(std::imag(d2psi_list[0].get()[1])) == 0.24972839540365);
692  CHECK(Approx(std::imag(dpsi_list[0].get()[1][0])) == -0.14785027769482);
693  CHECK(Approx(std::imag(dpsi_list[0].get()[1][1])) == 0.05546078377043);
694  CHECK(Approx(std::imag(dpsi_list[0].get()[1][2])) == -0.19276310593334);
695 
696  CHECK(Approx(std::imag(psi_list[1].get()[0])) == 0.00921878891696);
697  CHECK(Approx(std::imag(d2psi_list[1].get()[0])) == -0.01501853597019);
698  CHECK(Approx(std::imag(dpsi_list[1].get()[0][0])) == 0.00531194469490);
699  CHECK(Approx(std::imag(dpsi_list[1].get()[0][1])) == 0.00017539927704);
700  CHECK(Approx(std::imag(dpsi_list[1].get()[0][2])) == -0.00531168663225);
701  CHECK(Approx(std::imag(psi_list[1].get()[1])) == -0.06736883166719);
702  CHECK(Approx(std::imag(d2psi_list[1].get()[1])) == 0.19452200419216);
703  CHECK(Approx(std::imag(dpsi_list[1].get()[1][0])) == -0.13949226666533);
704  CHECK(Approx(std::imag(dpsi_list[1].get()[1][1])) == 0.05635953434634);
705  CHECK(Approx(std::imag(dpsi_list[1].get()[1][2])) == -0.17227552311686);
706  }
707 
708 
709  SECTION("compare MW/SW V/VGL")
710  {
711  //real
712  for (size_t iorb = 0; iorb < norb; iorb++)
713  {
714  CHECK(Approx(std::real(psi_list[0].get()[iorb])) == std::real(psiref_0[iorb]));
715  CHECK(Approx(std::real(psi_list[1].get()[iorb])) == std::real(psiref_1[iorb]));
716  CHECK(Approx(std::real(d2psi_list[0].get()[iorb])) == std::real(d2psiref_0[iorb]));
717  CHECK(Approx(std::real(d2psi_list[1].get()[iorb])) == std::real(d2psiref_1[iorb]));
718  for (size_t idim = 0; idim < SPOSet::DIM; idim++)
719  {
720  CHECK(Approx(std::real(dpsi_list[0].get()[iorb][idim])) == std::real(dpsiref_0[iorb][idim]));
721  CHECK(Approx(std::real(dpsi_list[1].get()[iorb][idim])) == std::real(dpsiref_1[iorb][idim]));
722  }
723  }
724  for (size_t iorb = 0; iorb < norb; iorb++)
725  for (size_t iw = 0; iw < nw; iw++)
726  CHECK(Approx(std::real(psi_v_list[iw].get()[iorb])) == std::real(psi_list[iw].get()[iorb]));
727  //imag
728  for (size_t iorb = 0; iorb < norb; iorb++)
729  {
730  CHECK(Approx(std::imag(psi_list[0].get()[iorb])) == std::imag(psiref_0[iorb]));
731  CHECK(Approx(std::imag(psi_list[1].get()[iorb])) == std::imag(psiref_1[iorb]));
732  CHECK(Approx(std::imag(d2psi_list[0].get()[iorb])) == std::imag(d2psiref_0[iorb]));
733  CHECK(Approx(std::imag(d2psi_list[1].get()[iorb])) == std::imag(d2psiref_1[iorb]));
734  for (size_t idim = 0; idim < SPOSet::DIM; idim++)
735  {
736  CHECK(Approx(std::imag(dpsi_list[0].get()[iorb][idim])) == std::imag(dpsiref_0[iorb][idim]));
737  CHECK(Approx(std::imag(dpsi_list[1].get()[iorb][idim])) == std::imag(dpsiref_1[iorb][idim]));
738  }
739  }
740  for (size_t iorb = 0; iorb < norb; iorb++)
741  for (size_t iw = 0; iw < nw; iw++)
742  CHECK(Approx(std::imag(psi_v_list[iw].get()[iorb])) == std::imag(psi_list[iw].get()[iorb]));
743  }
744  // app_log() << "vgl_refvalues: \n" << std::setprecision(14);
745  // for (int iorb = 0; iorb < 2; iorb++)
746  // {
747  // app_log() << "CHECK(Approx(std::real(psiref_0[" << iorb << "])) == " << std::real(psiref_0[iorb]) << ");\n";
748  // app_log() << "CHECK(Approx(std::real(d2psiref_0[" << iorb << "])) == " << std::real(d2psiref_0[iorb]) << ");\n";
749  // for (int idim = 0; idim < 3; idim++)
750  // app_log() << "CHECK(Approx(std::real(dpsiref_0[" << iorb << "][" << idim
751  // << "])) == " << std::real(dpsiref_0[iorb][idim]) << ");\n";
752  // }
753  // for (int iorb = 0; iorb < 2; iorb++)
754  // {
755  // app_log() << "CHECK(Approx(std::real(psiref_1[" << iorb << "])) == " << std::real(psiref_1[iorb]) << ");\n";
756  // app_log() << "CHECK(Approx(std::real(d2psiref_1[" << iorb << "])) == " << std::real(d2psiref_1[iorb]) << "));\n";
757  // for (int idim = 0; idim < 3; idim++)
758  // app_log() << "CHECK(Approx(std::real(dpsiref_1[" << iorb << "][" << idim
759  // << "])) == " << std::real(dpsiref_1[iorb][idim]) << ");\n";
760  // }
761  // app_log() << "vgl_refvalues: \n" << std::setprecision(14);
762  // for (int iorb = 0; iorb < 2; iorb++)
763  // {
764  // app_log() << "CHECK(Approx(std::imag(psiref_0[" << iorb << "])) == " << std::imag(psiref_0[iorb]) << ");\n";
765  // app_log() << "CHECK(Approx(std::imag(d2psiref_0[" << iorb << "])) == " << std::imag(d2psiref_0[iorb]) << ");\n";
766  // for (int idim = 0; idim < 3; idim++)
767  // app_log() << "CHECK(Approx(std::imag(dpsiref_0[" << iorb << "][" << idim
768  // << "])) == " << std::imag(dpsiref_0[iorb][idim]) << ");\n";
769  // }
770  // for (int iorb = 0; iorb < 2; iorb++)
771  // {
772  // app_log() << "CHECK(Approx(std::imag(psiref_1[" << iorb << "])) == " << std::imag(psiref_1[iorb]) << ");\n";
773  // app_log() << "CHECK(Approx(std::imag(d2psiref_1[" << iorb << "])) == " << std::imag(d2psiref_1[iorb]) << "));\n";
774  // for (int idim = 0; idim < 3; idim++)
775  // app_log() << "CHECK(Approx(std::imag(dpsiref_1[" << iorb << "][" << idim
776  // << "])) == " << std::imag(dpsiref_1[iorb][idim]) << ");" << std::endl;
777  // }
778  }
779 
780  SECTION("LCAOrbitalSet::mw_evaluateDetRatios")
781  {
782  // make VPs
783  const size_t nvp_ = 4;
784  const size_t nvp_2 = 3;
785  const std::vector<size_t> nvp_list = {nvp_, nvp_2};
786  VirtualParticleSet VP_(elec_, nvp_);
787  VirtualParticleSet VP_2(elec_2, nvp_2);
788 
789  // move VPs
790  std::vector<ParticleSet::SingleParticlePos> newpos_vp_(nvp_);
791  std::vector<ParticleSet::SingleParticlePos> newpos_vp_2(nvp_2);
792  for (int i = 0; i < nvp_; i++)
793  {
794  newpos_vp_[i][0] = 0.1 * i;
795  newpos_vp_[i][1] = 0.2 * i;
796  newpos_vp_[i][2] = 0.3 * i;
797  }
798  for (int i = 0; i < nvp_2; i++)
799  {
800  newpos_vp_2[i][0] = 0.2 * i;
801  newpos_vp_2[i][1] = 0.3 * i;
802  newpos_vp_2[i][2] = 0.4 * i;
803  }
804 
805  const auto& ei_table_ = elec_.getDistTableAB(ei_table_index);
806  const auto& ei_table_2 = elec_2.getDistTableAB(ei_table_index);
807  // make virtual move of elec 0, reference ion 1
808  NLPPJob<SPOSet::RealType> job_(1, 0, ei_table_.getDistances()[0][1], -ei_table_.getDisplacements()[0][1]);
809  // make virtual move of elec 1, reference ion 3
810  NLPPJob<SPOSet::RealType> job_2(3, 1, ei_table_2.getDistances()[1][3], -ei_table_2.getDisplacements()[1][3]);
811 
812  // make VP refvec
813  RefVectorWithLeader<VirtualParticleSet> vp_list(VP_, {VP_, VP_2});
814 
815  ResourceCollection vp_res("test_vp_res");
816  VP_.createResource(vp_res);
817  ResourceCollectionTeamLock<VirtualParticleSet> mw_vpset_lock(vp_res, vp_list);
818  VirtualParticleSet::mw_makeMoves(vp_list, p_list, {newpos_vp_, newpos_vp_2}, {job_, job_2}, false);
819 
820  // fill invrow with dummy data for each walker
821  OffloadVector<SPOSet::ValueType> psiMinv_data_0(norb), psiMinv_data_1(norb);
822  LCAOrbitalSet::ValueVector psiMinv_ref_0(psiMinv_data_0.data(), norb), psiMinv_ref_1(psiMinv_data_1.data(), norb);
823  for (int i = 0; i < norb; i++)
824  {
825  psiMinv_data_0[i] = 0.1 * i;
826  psiMinv_data_1[i] = 0.2 * i;
827  }
828  psiMinv_data_0.updateTo();
829  psiMinv_data_1.updateTo();
830 
831  std::vector<const SPOSet::ValueType*> invRow_ptr_list;
832  if (spo->isOMPoffload())
833  invRow_ptr_list = {psiMinv_data_0.device_data(), psiMinv_data_1.device_data()};
834  else
835  invRow_ptr_list = {psiMinv_data_0.data(), psiMinv_data_1.data()};
836 
837  // ratios_list
838  std::vector<std::vector<SPOSet::ValueType>> ratios_list(nw);
839  for (size_t iw = 0; iw < nw; iw++)
840  ratios_list[iw].resize(nvp_list[iw]);
841 
842  // just need dummy refvec with correct size
843  SPOSet::ValueVector tmp_psi0(norb), tmp_psi1(norb);
844  RefVector<SPOSet::ValueVector> tmp_psi_list{tmp_psi0, tmp_psi1};
845  spo->mw_evaluateDetRatios(spo_list, RefVectorWithLeader<const VirtualParticleSet>(VP_, {VP_, VP_2}), tmp_psi_list,
846  invRow_ptr_list, ratios_list);
847 
848  std::vector<SPOSet::ValueType> ratios_ref_0(nvp_);
849  std::vector<SPOSet::ValueType> ratios_ref_1(nvp_2);
850  // single-walker functions for reference
851  spo->evaluateDetRatios(VP_, tmp_psi0, psiMinv_ref_0, ratios_ref_0);
852  spo_2->evaluateDetRatios(VP_2, tmp_psi1, psiMinv_ref_1, ratios_ref_1);
853 
854  for (int ivp = 0; ivp < nvp_; ivp++)
855  CHECK(Approx(std::real(ratios_list[0][ivp])) == std::real(ratios_ref_0[ivp]));
856  for (int ivp = 0; ivp < nvp_2; ivp++)
857  CHECK(Approx(std::real(ratios_list[1][ivp])) == std::real(ratios_ref_1[ivp]));
858  for (int ivp = 0; ivp < nvp_; ivp++)
859  CHECK(Approx(std::imag(ratios_list[0][ivp])) == std::imag(ratios_ref_0[ivp]));
860  for (int ivp = 0; ivp < nvp_2; ivp++)
861  CHECK(Approx(std::imag(ratios_list[1][ivp])) == std::imag(ratios_ref_1[ivp]));
862 
863  CHECK(Approx(std::real(ratios_ref_0[0])) == 0.0963445284);
864  CHECK(Approx(std::real(ratios_ref_0[1])) == 0.0784621772);
865  CHECK(Approx(std::real(ratios_ref_0[2])) == 0.0312479567);
866  CHECK(Approx(std::real(ratios_ref_0[3])) == -0.0240189529);
867  CHECK(Approx(std::real(ratios_ref_1[0])) == 0.1926890568);
868  CHECK(Approx(std::real(ratios_ref_1[1])) == 0.1037508495);
869  CHECK(Approx(std::real(ratios_ref_1[2])) == -0.0747915097);
870 
871  CHECK(Approx(std::real(ratios_list[0][0])) == 0.0963445284);
872  CHECK(Approx(std::real(ratios_list[0][1])) == 0.0784621772);
873  CHECK(Approx(std::real(ratios_list[0][2])) == 0.0312479567);
874  CHECK(Approx(std::real(ratios_list[0][3])) == -0.0240189529);
875  CHECK(Approx(std::real(ratios_list[1][0])) == 0.1926890568);
876  CHECK(Approx(std::real(ratios_list[1][1])) == 0.1037508495);
877  CHECK(Approx(std::real(ratios_list[1][2])) == -0.0747915097);
878 
879  CHECK(Approx(std::imag(ratios_ref_0[0])) == -0.0090812301);
880  CHECK(Approx(std::imag(ratios_ref_0[1])) == -0.0385825385);
881  CHECK(Approx(std::imag(ratios_ref_0[2])) == -0.0610830209);
882  CHECK(Approx(std::imag(ratios_ref_0[3])) == -0.0809775403);
883  CHECK(Approx(std::imag(ratios_ref_1[0])) == -0.0181624602);
884  CHECK(Approx(std::imag(ratios_ref_1[1])) == -0.0856868673);
885  CHECK(Approx(std::imag(ratios_ref_1[2])) == -0.1487774316);
886 
887  CHECK(Approx(std::imag(ratios_list[0][0])) == -0.0090812301);
888  CHECK(Approx(std::imag(ratios_list[0][1])) == -0.0385825385);
889  CHECK(Approx(std::imag(ratios_list[0][2])) == -0.0610830209);
890  CHECK(Approx(std::imag(ratios_list[0][3])) == -0.0809775403);
891  CHECK(Approx(std::imag(ratios_list[1][0])) == -0.0181624602);
892  CHECK(Approx(std::imag(ratios_list[1][1])) == -0.0856868673);
893  CHECK(Approx(std::imag(ratios_list[1][2])) == -0.1487774316);
894 
895  // // print SW ref values
896  // for (int ivp = 0; ivp < nvp_; ivp++)
897  // app_log() << "CHECK(Approx(std::real(ratios_ref_0[" << ivp << "])) == " << std::real(ratios_ref_0[ivp]) << ");\n";
898  // for (int ivp = 0; ivp < nvp_2; ivp++)
899  // app_log() << "CHECK(Approx(std::real(ratios_ref_1[" << ivp << "])) == " << std::real(ratios_ref_1[ivp]) << ");\n";
900 
901  // // print MW ref values
902  // for (int iw = 0; iw < nw; iw++)
903  // for (int ivp = 0; ivp < nvp_list[iw]; ivp++)
904  // app_log() << "CHECK(Approx(std::real(ratios_list[" << iw << "][" << ivp
905  // << "])) == " << std::real(ratios_list[iw][ivp]) << ");\n";
906  // for (int iw = 0; iw < nw; iw++)
907  // for (int ivp = 0; ivp < nvp_list[iw]; ivp++)
908  // app_log() << "CHECK(Approx(std::imag(ratios_list[" << iw << "][" << ivp
909  // << "])) == " << std::imag(ratios_list[iw][ivp]) << ");\n";
910  }
911 }
912 
913 TEST_CASE("LCAOrbitalSet batched PBC DiamondC", "[wavefunction]")
914 {
915  SECTION("2x1x1 real") { test_LCAO_DiamondC_2x1x1_real(false); }
916 #ifdef QMC_COMPLEX
917  SECTION("2x1x1 cplx") { test_LCAO_DiamondC_2x1x1_cplx(false); }
918 #endif
919 }
920 
921 TEST_CASE("LCAOrbitalSet batched PBC DiamondC offload", "[wavefunction]")
922 {
923  SECTION("2x1x1 real offload") { test_LCAO_DiamondC_2x1x1_real(true); }
924 #ifdef QMC_COMPLEX
925  SECTION("2x1x1 cplx offload") { test_LCAO_DiamondC_2x1x1_cplx(true); }
926 #endif
927 }
928 } // namespace qmcplusplus
qmcplusplus::CrystalLattice
a class that defines a supercell in D-dimensional Euclean space.
Definition: CrystalLattice.h:55
qmcplusplus::ParticleSet::setName
void setName(const std::string &aname)
Definition: ParticleSet.h:237
Libxml2Document
class that handles xmlDoc
Definition: Libxml2Doc.h:76
qmcplusplus::SpeciesSet::addSpecies
int addSpecies(const std::string &aname)
When a name species does not exist, add a new species.
Definition: SpeciesSet.cpp:33
qmcplusplus::Units::real
QMCTraits::RealType real
Definition: unit_conversion.h:23
qmcplusplus
helper functions for EinsplineSetBuilder
Definition: Configuration.h:43
qmcplusplus::VirtualParticleSet::createResource
void createResource(ResourceCollection &collection) const
initialize a shared resource and hand it to a collection
Definition: VirtualParticleSet.cpp:76
qmcplusplus::app_log
std::ostream & app_log()
Definition: OutputManager.h:65
qmcplusplus::LatticeParser::put
bool put(xmlNodePtr cur)
Definition: LatticeIO.cpp:32
qmcplusplus::TEST_CASE
TEST_CASE("complex_helper", "[type_traits]")
Definition: test_complex_helper.cpp:38
qmcplusplus::VirtualParticleSet
A ParticleSet that handles virtual moves of a selected particle of a given physical ParticleSet Virtu...
Definition: VirtualParticleSet.h:39
qmcplusplus::Vector::device_data
pointer device_data()
Return the device_ptr matching X if this is a vector attached or owning dual space memory...
Definition: OhmmsVector.h:245
Libxml2Document::getRoot
xmlNodePtr getRoot()
Definition: Libxml2Doc.h:88
qmcplusplus::ParticleSet::createResource
void createResource(ResourceCollection &collection) const
initialize a shared resource and hand it to a collection
Definition: ParticleSet.cpp:944
qmcplusplus::test_LCAO_DiamondC_2x1x1_cplx
void test_LCAO_DiamondC_2x1x1_cplx(const bool useOffload)
Definition: test_LCAO_diamondC_2x1x1.cpp:454
qmcplusplus::LCAOrbitalBuilder
SPOSetBuilder using new LCAOrbitalSet and Soa versions.
Definition: LCAOrbitalBuilder.h:33
qmcplusplus::ParticleSet::update
void update(bool skipSK=false)
update the internal data
Definition: ParticleSet.cpp:349
OHMMS::Controller
Communicate * Controller
Global Communicator for a process.
Definition: Communicate.cpp:35
Array::resize
void resize(const std::array< SIZET, D > &dims)
Resize the container.
Definition: OhmmsArray.h:65
qmcplusplus::SpeciesSet::addAttribute
int addAttribute(const std::string &aname)
for a new attribute, allocate the data, !More often used to get the index of a species ...
Definition: SpeciesSet.cpp:45
qmcplusplus::ParticleSet::getDistTableAB
const DistanceTableAB & getDistTableAB(int table_ID) const
get a distance table by table_ID and dyanmic_cast to DistanceTableAB
Definition: ParticleSet.cpp:344
Communicate
Wrapping information on parallelism.
Definition: Communicate.h:68
qmcplusplus::lattice
CrystalLattice< OHMMS_PRECISION, OHMMS_DIM > lattice
Definition: test_SpinDensityNew.cpp:142
qmcplusplus::ParticleSet::addTable
int addTable(const ParticleSet &psrc, DTModes modes=DTModes::ALL_OFF)
add a distance table
Definition: ParticleSet.cpp:307
qmcplusplus::ParticleSet
Specialized paritlce class for atomistic simulations.
Definition: ParticleSet.h:55
qmcplusplus::imag
float imag(const float &c)
imaginary part of a scalar. Cannot be replaced by std::imag due to AFQMC specific needs...
Definition: complex_help.hpp:91
qmcplusplus::test_LCAO_DiamondC_2x1x1_real
void test_LCAO_DiamondC_2x1x1_real(const bool useOffload)
Definition: test_LCAO_diamondC_2x1x1.cpp:37
qmcplusplus::QMCTraits::ValueType
QTBase::ValueType ValueType
Definition: Configuration.h:60
qmcplusplus::REQUIRE
REQUIRE(std::filesystem::exists(filename))
qmcplusplus::SPOSet::ValueVector
OrbitalSetTraits< ValueType >::ValueVector ValueVector
Definition: SPOSet.h:49
qmcplusplus::LCAOrbitalBuilder::createSPOSetFromXML
std::unique_ptr< SPOSet > createSPOSetFromXML(xmlNodePtr cur) override
create an sposet from xml (legacy)
Definition: LCAOrbitalBuilder.cpp:471
qmcplusplus::ParticleSet::R
ParticlePos R
Position.
Definition: ParticleSet.h:79
qmcplusplus::VirtualParticleSet::mw_makeMoves
static void mw_makeMoves(const RefVectorWithLeader< VirtualParticleSet > &vp_list, const RefVectorWithLeader< ParticleSet > &p_list, const RefVector< const std::vector< PosType >> &deltaV_list, const RefVector< const NLPPJob< RealType >> &joblist, bool sphere)
Definition: VirtualParticleSet.cpp:177
qmcplusplus::ParticleSet::getSpeciesSet
SpeciesSet & getSpeciesSet()
retrun the SpeciesSet of this particle set
Definition: ParticleSet.h:231
qmcplusplus::ParticleSet::create
void create(const std::vector< int > &agroup)
create grouped particles
Definition: ParticleSet.cpp:125
qmcplusplus::RefVector
std::vector< std::reference_wrapper< T > > RefVector
Definition: template_types.hpp:32
qmcplusplus::ParticleSet::print
void print(std::ostream &os, const size_t maxParticlesToPrint=0) const
print particle coordinates to a std::ostream
Definition: ParticleSet.cpp:275
qmcplusplus::SPOSet::GradVector
OrbitalSetTraits< ValueType >::GradVector GradVector
Definition: SPOSet.h:51
Libxml2Document::parseFromString
bool parseFromString(const std::string_view data)
Definition: Libxml2Doc.cpp:204
qmcplusplus::ParticleSet::mw_update
static void mw_update(const RefVectorWithLeader< ParticleSet > &p_list, bool skipSK=false)
batched version of update
Definition: ParticleSet.cpp:362
qmcplusplus::CHECK
CHECK(log_values[0]==ComplexApprox(std::complex< double >{ 5.603777579195571, -6.1586603331188225 }))
qmcplusplus::ResourceCollectionTeamLock
handles acquire/release resource by the consumer (RefVectorWithLeader type).
Definition: ResourceCollection.h:58
qmcplusplus::NLPPJob
meta data for NLPP calculation of a pair of ion and electron This is not just meta data...
Definition: VirtualParticleSet.h:30
qmcplusplus::LCAOrbitalSet
class to handle linear combinations of basis orbitals used to evaluate the Dirac determinants.
Definition: LCAOrbitalSet.h:30
WaveFunctionComponent.h
Declaration of WaveFunctionComponent.
qmcplusplus::SpeciesSet
Custom container for set of attributes for a set of species.
Definition: SpeciesSet.h:33
qmcplusplus::Vector::updateTo
void updateTo(size_type size=0, std::ptrdiff_t offset=0)
Definition: OhmmsVector.h:263
Array
A D-dimensional Array class based on PETE.
Definition: OhmmsArray.h:25